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ABSTRACT 


Significant problems are presented in the vibration and rotation 
analysis of spacecraft with distributed structural flexibility and 
momentum exchange controllers. These systems exhibit gyroscopic cou- 
pllng which depends on the rotor speed and orientation, which must 
remain explicit in the analysis as control variables. 

An investigation is made into floating reference frames in order 
to allow first order vibration analysis in the presence of large system 
rotations. When the deformations of an elastic continuum are expanded 
in terms of the free-free modes of an unconstrained system, the rigid 
body modes are found to be fixed relative to the Tisserand frame, with 
respect to which the relative momentum is zero. The proof presented 
for this is based on the orthogonality condition for modes with distinct 
natural frequencies. This result also guarantees the Independence of 
coordinates for all modes with nonzero natural frequencies. A Modified 
Tisserand Constraint is introduced in order to define a floating refer- 
ence frame with similar properties for an elastic body which contains a 
spinning rotor. 

Finite element equations of motion are derived for a completely 
flexible spacecraft with momentum exchange controllers, using a Modified 
Tisserand Frame. The deformable systems covered in this application 
are assumed to undergo only • small rotations, and therefore the rotor 
torques must formally be small, although in engineering applications it 
may be possible to relax this constraint. A modal analysis is performed 
for the system and the resulting set of equations is reduced in number 
by a truncation procedure for more efficient system simulation. 
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In order to gain Insight Into system behavior, a continuum 
analysis is performed for a simple physical system consisting of a 
uniform beam and an axisymmetric rotor,. Equations are derived using 
Hamilton's principle and closed form solutions are obtained using 
generalized methods of separation of variables. 

Numerical examples are presented for both the finite element and 
continuxjm equations of motion. 
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CHAPTER 1 


INTRODUCTION 

1.1 Statement of Problem 

Modern spacecraft design relies on the successful dynamic analysis 
of the system In order to provide attitude stabilization. This Is true 
for both active and passive attitude control systems. Present systems 
emphasize active control techniques with momentum exchange controllers 
being the most popular because of their greater accuracy potential. 

Momentum exchange controllers fall Into three general categories; 
dual spin, momentum (or reaction) wheel and control moment gyro systems. 
The dual spin spacecraft typically has a large rotor spinning at a 
nearly constant speed and a despun platform. For communications satel- 
lites, the antennas are mounted on the despun platform, which Is aimed 
at a point on the earth. The platform orientation about one axis Is 
controlled by changes In the rotor speed, and orientation about the 
other two axes Is passively controlled by a nutation damper. An example 
of this type of system Is the Intelsat IV spacecraft. The second type 
Is the momentum wheel or reaction wheel system. The reaction wheel has 
a nominal spin rate near zero. There may be three or more rotors, and 
control about each axis Is accomplished actively by controlling rotor 
speeds. This Is often referred to as a type of a three axis stabilized 
system. An example of this system is the Orbiting Astronomical Obser- 
vatory (OAO) . The third type of system uses the control moment gyro 
(CMC). This device consists of a glmballed constant speed rotor. 

Torques are applied to the vehicle by driving the gimbals and changing 
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the orientation of the rotor. An example of a spacecraft with this 
system is SKYLAB. A mixture of the above and active control using gas 
jets has also been Incorporated in spacecraft design. 

To analyze these systems, many powerful techniques and formulations 
are available. If the system is treated as a group of interconnected 
rigid bodies, Newton's laws may be used to solve for the equations of 
motion. Advanced formulations exist for general classes of configura- 
tions which use nested body techniques to eliminate the constraint 
forces. An example of this is the Hooker-Margulles-Hooker system of 
equations. If the bodies are not rigid, flexibility of terminal append- 
ages may also be analyzed by using the "hybrid coordinate" method. This 
approach utilizes a combination of discrete and modal coordinates. The 
underlying finite element formulation allows truncation of the number of 
degrees of freedom after a modal analysis and therefore has the very 
great advantage of computational ease. Appendix D contains a simpli- 
fied example of this approach. 

Ongoing developments are causing the rigid body model of a space- 
craft to be less acceptable. The effects of distributed structural 
flexibility need to be included to reduce spacecraft weights and to meet 
ever more stringent accuracy requirements. The increase in size of 
spacecraft due to the advent of the Space Shuttle will also make rigid 
body models less precise. This research was motivated by the attitude 
control problems associated with the proposed Large Space Telescope 
(LST) . 

The incorporation of distributed structural flexibility for the 
entire vehicle introduces great complexity to the dynamic analysis. The 


linearized equations of motion for a finite element model with momentum 
exchange controllers contain time varying coefficients. Thus, an eigen- 
value analysis to find vehicle normal modes is not possible. If the 
spin speeds are constant, a modal analysis results in complex modes; 
that is, the eigenvectors are represented by complex numbers. If a 
continuum model is adopted, the resulting partial differential equations 
are coupled and have time varying boundary conditions. 

1. 2 Motivation and Objectives 

The motivation for this study was a desire to address the 
difficulties posed by spacecraft of the next level of complexity. A 
review of the present techniques does not indicate an appropriate 
course of action for the dynamic analysis of spacecraft with both dis- 
tributed flexibility and momentum exchange controllers. A finite ele- 
ment model must be adopted to allow for arbitrary mass and stiffness 
properties. A modal analysis will allow the introduction of distributed 
coordinates and the truncation to a smaller number of coordinates. A 
modal analysis can be accomplished only for nearly constant speeds, and 
then the complex mode shapes depend on the chosen spin speed and rotor 
orientation. A method of Staging mode shapes for different ranges of 
spin speeds or orientations would have to be devised to handle a general 
system. This is not an attractive method. 

The hybrid coordinate method is frustrated by the lack of a 
central rigid body. For this method, the modal coordinates are coupled 
to each other through the motion of the central rigid body. If a rigid 
body is designated from among the rigid bodies in the finite element 
model, perhaps the body to which the rotor is mounted, then the mass 
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properties of the body will Influence the determination of the modal 
coordinates. This casts the validity of the truncation procedure into 
doubt . 

The principal objectives of this research on the vibration of 
spacecraft with distributed flexibility and momentum exchange control- 
lers are: 

1. Provide the advantage of a first order vibration analysis 

by adopting a floating reference frame which somehow "moves" 
with the flexible body. 

2. Develop equations of motion for use in spacecraft dynamics 
and control analysis. 

3. Gain insight into system behavior through analysis of simple, 
physical systems. 

1.3 Scope of the Dissertation 

The body of the work pursued in the dissertation is broken into 
two main parts. The first part consists of the introduction and 
analysis of various floating reference frames. The second part consists 
of the formulation of equations of motion for the class of spacecraft 
characterized by distributed flexibility and momentum exchange control- 
lers, using both finite element and continuum models, numerical 
examples are developed for both classes of models, using eigenvalue 
analysis and numerical integration for the finite element model and 
closed form solutions for the continuum model . The two parts of the 
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dissertation are interrelated by the use of a floating reference frame 
in the finite element analysis of a system with momentum exchange con- 
trollers. In order to aid in understanding the organization of the 
dissertation, a brief summary of each chapter is presented. 
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Chapter one introduces the problem area and presents the 
motivation for the studies undertaken. The difficult aspects of the 
problem are identified and the existing solution techniques are 
reviewed. The topics covered by the dissertation are discussed. 

Chapter two introduces the concept of a floating reference frame. 
For an elastic body free to rotate, these frames are Introduced in 
order to ensure the validity of a first order vibration analysis. A 
review is made of several distinct ways of defining a floating reference 
frame. They all provide for frame motion which somehow "follows" the 
overall body. In order to demonstrate how this is accomplished, a very 
simple example problem is worked out to Illustrate the motion of 
specific frames when the body is deformed. The chapter concludes with 
a brief historical perspective of the use of floating frames to solve 
problems in the physical sciences. 

Chapter three deals with a specific floating frame, the Tisserand 
frame, which is defined in terms of angular momentum relative to the 
frame. An examination of the constraint relationship reveals a rather 
elegant solution form which Involves the free-free modes of an uncon- 
strained system. This is called the Mode Shape Constraint. An example 
is worked out in order to Illustrate the application of the Mode Shape 
Constraint to a simple problem. The last section of the chapter expands 
the concept of the Tisserand frame to systems containing spinning rigid 
rotors. This is done with the introduction of the Modified Tisserand 
Constraint. 

Chapter four begins the second main part of the dissertation. 

The focus now turns to the derivation of equations of motion for systems 
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containing a spinning rigid rotor. System behavior is limited by the 
requirement that variations in rotor angular momentum be first order. 

The deformable body is also assumed to undergo only small, first order 
rotations. This chapter uses a Finite Element Model and applies the 
Modified Tisserand Constraint. A modal analysis leads to a system of 
equations which are still coupled, but which may be truncated to a 
smaller number of coordinates. 

Chapter five defines a very simple physical system consisting of 
a uniform beam and axisymmetric rotor for analysis as a continuum. 
Hamilton's Principle is used to derive the governing partial differen- 
tial equations with boundary conditions. Generalized separation of 
variables techniques are used in order to solve for the system natural 
frequencies and system normal modes. Numerical solutions are provided 
for a specific choice of parameters in order to illustrate the system 
behavior. 

Chapter six presents numerical results for the Finite Element 
Model equations of motion. The strategy is quite simple. The simple 
physical system Introduced in Chapter five is modeled by finite elements 
and initial conditions are chosen for a vehicle normal mode as deter- 
mined by the continuum analysis. Numerical integration of the Finite 
Element Model equations of motion should then demonstrate periodic 
behavior corresponding to the system natural frequency. The method of 
truncation to a smaller number of modal coordinates is investigated. 

Chapter seven summarizes the dissertation problem, techniques and 
contributions. The principal conclusions of the research are assembled 
and areas of future work are outlined . 
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CHAPTER 2 


FLOATING REFERENCE FRAMES 

2.1 Chapter Summary 

The problem of specifying a set of axes for a deformable body is 
a difficult one. The use of a moving or ’’floating" set of axes is 
introduced here. The choice of a floating reference system is not 
unique and several alternative methods are described. In particular, 
there are the locally attached frame, the principal axis frame, the 
Tisserand and Buckens frames, and the rigid body mode frame. The 
advantages of each different type of floating frame are discussed and 
a simple example is used to illustrate the motion of the different 
floating frames. The historical origin of the principal axis and 
Tisserand reference frames is briefly discussed to demonstrate their 
application to past problems in the physical sciences. 

2.2 Floating Frames for Deformable Bodies 

For the problems discussed here, we are concerned with small 
deformations of elastic bodies. Specifically, we will define a deform- 
able body to be a body for which the relative displacements are so 
small that only first order terms need be retained in the analysis. 

The body is allowed overall motion which is completely unrestricted; 
that is, the body has complete freedom of motion in responding to 
impressed moments and forces. 

For such a body, it is difficult to specify a set of axes from 
which to measure deformations. If an Inertially fixed set of axes is 
chosen, the displacements relative to these axes may grow large if the 
body undergoes any appreciable rotation due to an externally applied 
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. Then the dynamic analysis of the system using these 
displacements would require more than a first order analysis. In order 
to simplify the analysis, the idea presents itself to use a frame that 
somehow moves with the body. If a frame moves with the body, or 
"floats", in the proper way, then the displacements measured relative 
to this floating frame will be small. The dynamic analysis of the 
system may then be pursued using a first order analysis. It is this 
property, the use of a first order analysis, that makes a floating 
frame attractive. 

There are five types of floating frames that will be treated 
here; they are: 

• Locally attached frame 

f Principal axis frame 

• Tisserand frame 

• I Buckens frame 

• Rigid body mode frame. 

Before these different types of frames are discussed in detail, several 
general characteristics of floating frames will be discussed. A 
deformable body experiences relative displacement and therefore the 
system inertia quantities do hot remain constant. In the formulation 
of equations of motion, the system angular momentum and kinetic energy 
are a function of the relative displacements and possess a more com- 
plicated structure than for rigid bodies. The system angular momentum 
for a deformable body about its center of mass is written as; 

Ji - / P ^ P dm (2.1) 

•'b 
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where 


£ = vector from the center of mass to a generic mass element " 

• _ 

£ " dt £ “ derivative of £ relative to Inertial space 
dm = generic mass element 

D = integration over the deformable body 

H =* angular momentum of the system about the center of mass. 

If a floating frame, f, is Introduced with its origin at the center of 
mass, then 

i‘i+““xp (2.2) 

where 


£ “ dt £ “ derivative of £ relative to floating frame f . 

,,fi , ^ 

0) = angular velocity of frame f relative to inertial space. 

Using this relationship to evaluate the angular momentum yields 

H = [] • + ^£ X P dm (2.3) 

since 

. . ! i i ■ 

jTpx dn= Q- 0 )“ (2.4) 

where Q is the inertia dyadic for the mass center. The second quan- 
tity in Equation (2.3) is the angular momentum relative to the floating 
frame. This is referred to as the internal angular momentum. The 
first term of Equation (2.3) is structurally identical to the rigid 
body angular momentum, but since the body is deformable, the Inertia 
dyadic is not constant. It should be noted that for small (first order) 
displacements relative to the floating frame, the variations in the 
inertia dyadic will be first and higher order terms. 
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The kinetic energy of a body with an inertially fixed center of 
mass may be written 

T = “ / £ • £ dm . (2.5) 

n) 

Using the relationship for the time derivative in Equation (2.2) allows 
the kinetic energy to be rewritten using a floating frame as 

T - i • n • • f 

^ •'d 

where 

• n • f ^ p) • (£^^ ^ dm. 

^ : *D 

The first term of this equation has the same structural form as the 
kinetic energy of a rigid body. Again, this quantity differs from the 
rigid body quantity because the inertia dyadic is not constant. The 
second term includes the same internal angular momentum expression 
found earlier in Equation (2.3). The last term could be called the 
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internal kinetic energy, since it represents kinetic energy contributed 
by the dot product of velocities relative to the frame. 

2.3 Description of Specific Floating Frames 

This next section will now deal with the specific definitions for 
the most commonly available floating frames. The first type to be 
discussed is the locally attached frame. For this frame, a sub-body 
or mass element is identified in the deformable body, and a frame is 

defined that follows the motion of this sub-body or mass element. An 

; i ^ 

example of this type of frame may be associated with a spacecraft with 
a rigid central body and a flexible appendage. A reference frame may 
be attached to the central rigid body. Also, if the appendage is 
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(£ X £) dm + ^ J*£_ • £ dm (2.6) 


driven, that is, its orientation relative to the central rigid body 
varies as a specific function of time in addition to small displace- 
ments, a floating reference frame may be attached to a mass element in 
the appendage. This would most commonly be a set of axes fixed in the 
appendage at the mounting point. The angular momentum and kinetic 
energy for a locally attached frame would be given by equations (2.3) 
and (2.6) respectively. As a general rule, no simplification of these 
expressions could be guaranteed for the locally attached reference 
frame. 

The next type of floating reference, the principal axis frame, 
does offer some simplification of the expressions for angular momentum 
and kinetic energy. For this frame, the origin of the axes would be 
the center of mass and the orientation of the axes would be such that 
the Inertia matrix of the deformable body would be diagonal. This 
defines the location of the principal axes of inertia. The components 
of a diagonal inertia matrix are called the moments of inertia for 
principal axes and the products of inertia are all zero. Because the 
Inertia matrix contains only three components, the calculations for 
angular momentum and kinetic energy are simplified. However, this is 
done at the expense of introducing three constraint relationships that 
require the products of inertia all to be zero. The constraint rela- 
tionships are 


12 “ ^21 " “ J 

I— ! 

Q. 

\ 

H 

dm = 0 


13 ^31 ^ " J 


P 3 

dm = 0 


23 " ^32 " " J 


P 3 

dm = 0 

(2.7) 
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where (p^, ^ 2 * P 3 ^ components of £ in the floating 

frame 

For the Tisserand frame, the expressions for angular momentum 
and kinetic energy are structurally simplified by moving the axes so 

i '■ ^ ! i 

as to set the internal angular momentum always to zero. The require- 
ment is also made that the internal linear momentum be zero. This is 
accomplished by the simple requirement that the origin of the frame 
be located at the center of mass. In order to set the internal angular 
momentum^to zero, a constraint relationship is introduced: 

J £ X £ dm = 0 (2.8) 

•'d 

Referring back to Equation (2.3), the angular momentum is thus 

H = [] • (2.9) 

This is structurally identical to the rigid body form, although as 
noted, the inertia matrix is not a constant for a deformable body as it 
would be for a rigid body relative to axes fixed in that body. Going 
back to the constraint relationship, let us evaluate this in greater 
detail. First we introduce a new expression for the position relative 
to the center of mass 

£ = £ + £ ( 2 . 10 ) 
where £ is the position of a generic mass element in the undeformed 
state. This vector, £, is fixed in the floating reference frame and 
may be thought of as the Station location of a mass element. The 
vector u. represents deformation of a generic mass element. For a 
deformable body, this will be a first order quantity. The derivative 
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of the position of a generic mass element relative to the center of 
mass is then written as 

O * 0 o 

£=£+U=U (2.11) 

since £ is fixed in the frame and has no derivative relative to the 

frame. The constraint relationship is then 

J | £ X u dm = 0 . 

D 

This is equivalent to the three scalar equations 


( 2 . 12 ) 


J\ ^2^ ~ ^2^^ ^ dm = 0 


u^(z + u^) - u^Cx + u^) dm - 0 

J' u^Cx + u^) - u^(y + U 2 > dm = 0 


( 2 . 13 ) 


where (x, y, z) are components of £ in the floating frame. The next 
step in our examination of the constraint relationship is to Introduce 
a separation of variables for the deformations. Then the deformations 
may be written as 


n 


“i = X n. (t) 

i»l ^ ^ 


n . _ 

^2 “ X 4>o(p) n>(-t) 
i»i ^ ^ 


n 


^3 = X ‘i’o® Tli(t) 
i=l 


( 2 . 14 ) 


where the variables ((j)^, c})^, (J)^) are the components of the mode shape ^ 
and depend only on the spatial variables. The variable is the modal 
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coordinate and depends on the time. This separation of variables may 
be written In vector form as 

.1. 


ii - E in. 
1-1 ^ 


(2.15) 

Using the separation of variables In the constraint relationship yields 




dm rij = 0 


(2.16) 


This is a vector constraint of the Pfaffian form 


n 


2 Cn) n = 0 

j=l ^ 


(2.17) 


where 


Sj(n) 




(j)^ dm. 


(2.18) 


If the constraint relationship (Equation (2.16)) possesses an integral* 
then the constraint is holonomic. It is important to note that if a 
Tisserand constraint is used in a Lagrangian formulation, then the 
equations of motion must include Lagrange multipliers for the general 
case. The Tisserand frame is a general concept and may be applied to a 
system where the deformations are large. In the case of a deformable 
body as previously defined, the relative displacements are small and 
may be treated analytically as first order quantities. The analysis 
may then proceed using first order quantities and may ignore second and 
higher order terms. It is proceeding in this direction which leads us 
to introduce the Buckens frame. The Buckens constraint relationship 
is simply a first order Tisserand constraint. It may be written as 



( 2 . 19 ) 


^ £ X u dm s= 0 , 

Again the origin of the system is placed at the center of mass. The 
Buckens constraint is identical to the Tisserand constraint if this 
relationship holds 

• , J 

_u X ju dm = 0 (2.20) 

The justification for the above is that second order quantities may be 
ignored for a deformable body. Introducing a separation of variables 
again yields a constraint of the Phafflan form 
n 

T. 

j'l 

The difference is that now the coefficients of the constraint relation- 
ship are constant and the expression always possesses an integral. 

Thus, the Buckens or first order Tisserand constraint is a holonomic 
constraint. Since the Buckens and Tisserand constraints are identical 
for a deformable body as defined above, they will both be referred to 
as the Tisserand constraint. This choice is made because of the 
historical precedence of the Tisserand constraint and the fact that it 
is more common in the literature. This is not to minimize the contri- 
bution made by Buckens, for it is the first order form which will allow 
the greater use of the Tisserand frame that is explored in the next 
chapter. 

The last frame to be described is the rigid body mode frame. 

This concept arises in structural dynamics for semidefinite systems. 

A semidefinite system is one for which the strain energy may be zero 


J. 


p X (j)' 


dm T|j = 0 


( 2 . 21 ) 
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without the motion being zero. -Unrestrained systems, or systems 
without supports, are typical examples of semldefinite systems. A 
rigid body mode is defined as a displacement which results In zero 
strain energy. The rigid body mode frame follows the displacement 
which results in zero strain energy. There is an essential difference 
between the rigid body mode and the equilibrium position of a body. 

The rigid body mode is associated with zero strain energy which results 
in the equilibrium position of the body at rest. The body must be 
inertlally at rest both translationally and rotationally. If a body 
is spinning at a given rate, the equilibrium position will not coincide 
with the rigid body mode (equilibrium at rest) since the centrifugal 
forces may be thought of as inducing a nonzero strain energy at equilib- 
rium. The presence of environmental forces such as thermal gradients 
can complicate the notion of a rigid body mode. The dictum of zero 
strain energy must still be followed in this case, but the effects of 
environmental forces on the strain energy can be neglected. 

The constraint relationship associated with the rigid body mode 
is simply that the strain energy be zero. The difficulty of working 
with this constraint is circumvented by showing a relationship between 
the Tlsserand frame for deformable bodies and the rigid body mode frame. 
This allows the rigid body mode to be easily applied. This is pursued 
in the next chapter. 

2.4 Example of Frame Motion 

In order to give insight into the concept of a floating reference 
frame, a simple example problem has been formulated. This example sys- 
tem is shown in Figure 2.1 and consists of two uniform rigid bodies 
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that are connected by a line hinge. All motion is planar. In the 
undeforiiied position, the two bodies are aligned in a straight line. 

The length of body one is £ and the length of body two is The 

mass per unit length of body two is twice that of body one. This 
places the center of mass for the undeformed body at the hinge point. 
The above system is made non-uniform in order to separate the motion 
of the principal axis and Tisserand frames. If the bodies were iden- 
tical, the two frames would move in the same manner. The coordinate 
systems for the example are shown in Figure 2.2. The deformation of 
the system is specified by the angle a between the two bodies. A set 

/*v 'N ^ 

of axes, with respect to which unit vectors b^, and b^ are fixed, 
is fixed in body two with its origin at the hinge point. The location 
of the hinge point relative to the center of mass is specified by the 

/S /Vi 

vector R. A floating frame, in which unit vectors _f 2 , and 
fixed, has its origin at the center of mass and makes an angle 3 with 
the body fixed axes. A mass element is located by a vector x relative 
to the hinge point with the subscript specifying the associated body. 
The center of mass is defined by 

P dm = 0 (2.22) 



This becomes for the example 



(R 


4- x^) m dx + 


Jt 


■vT 


(jR + X2^ 



(2.23) 


where 

/s. 

- x(cosa + sina ^ 2 ) 

(2.?A) 
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This yields the position vector of the hinge point relative to the 
center of mass 


R = ]b^j^ (1 - cosa^ - h 2 ^-|- sinotj 


(2.25) 

Note that the center of mass is at the hinge point for zero deforma- 
tions . 

The Tisserand frame will have its origin at the center of mass 
and will obey the constraint relationship 


|pxp 


dm = 0 


(2.26) 


where 


£l R 

£-2 “ ~ ~2 * 


The constraint then becomes 


(2.27) 


1 


I “(Vly ■%%)'== ^/o ( 2 . 28 ) 

Evaluating the vectors in Equation (2.27) in the floating reference 
frame so as to facilitate the derivatives involved yields 


£l = W' 


^ Q ^ 

cosp + 


- Sin3 + 


(==-1) 
(=< - 1) 


cos (ot-3) 
ln(a-3) 


£2 - 


-^x - -|•jcds3 - cos(a-3) 

- ■^jsin3 ^ sin(a-3) 


(2.29) 


where {_f} is the 3 x 1 array of vectors f^^. 
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With these relationships and the derivatives of the components, the 
constraint relationship may be evaluated. After the necessary algebra 
and integrations. Equation (2,28) becomes 

(.2342 + .0991 cosa)a - (.3333 + .1982 cosa)3 =0 (2.30) 

This constraint is in the Pfaffian form and since it possesses an 
Integral, the Tisserand constraint here is holonoraic. If the deforma- 
tions are small, then the constraint shows that the frame must move 
relative to the bodies in a specific manner, that is 

3 = .586 a (2.31) 

For a Buckens frame, the constraint relationship is 


£ X dm = 0 

where 





The constraint can be stated in the form 



(2.32) 


(2.33) 


Jl 

- ■ 

I ra Pi p -1 dx +1 2m p^ P 2 dx = 0 (2.34) 

Jq -^x -^y Jq X y 

Evaluating the derivatives from Equation (2.29) and performing the 
necessary algebra gives the constraint in the form 

cos(a-3)(d-3) - .7071 cos33 = 0 (2.35) 

For small deformationp, the relationship between the angles will be 

3 = .586 a (2.36) 


This is the same answer for small deformations that resulted from the 
Tisserand frame. This is to be expected since the Buckens frame is 
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simply a first-order Tisserand frame and they are identical for systems 
undergoing small deformations . 

For the principal axis frame, the constraint will set the product 
of Inertia to zero 


■12 


■I 


p p dm - 0 

X y 


(2.37) 


The constraint relationship for the example is then written 

Jl ' 

pi \/2 

I mp,p, dx+l 2mp«p„ dx=0 

Jq -^y Jo 


(2.38) 


Using the expressions found in Equation (2.29) which hold for any float- 
ing frame, and completing the necessary algebra yields 
-.1367 sin3 cos$ -I- .0991 sin(a-23) + .2341 sin(a-3) cos(a-3) = 0 (2.39) 
For small angles, the relationship between the coordinates is 

3 = .476 a (2.40) 


The above constraint relationships are holonomlc in form since they 
relate the coordinates. The Tisserand and Buckens constraints are of 
the Pfafflan form since they relate velocities. They may be integrated 
in this case and if both angles are zero initially, then they yield 


3 = .586 a 


(2.41) 


This means that the principal axis frame is definitely distinct from 
the Tisserand frame for this example.’ The Tisserand frame will move 
more relative to the body fixed frame due to a deformation than w^ll 
the principal axis frame for this specific case. The Tisserand frame 
will rotate 1.23 degrees for each degree that the principal axis frame 
rotates due to a deformation of the system. The rotation of the frames 
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discussed above is relative to the body fixed frame. One special case 
of the rotation of the Tisserand frame relative to inertial space 
should be mentioned. For a moment free body with zero net angular 
velocity, the angular momentum will remain zero due to conservation 
laws. For this case, by reference to Equation (2.9), the frame angular 
velocity must be zero since the inertia matrix is positive definite. 
Thus, the Tisserand frame in this special case will be inertially fixed 
while the principal axis frame moves. 

2.5 Historical Perspective 

The principal axis and Tisserand frames both had applications to 
problems in the physical sciences during the late nineteenth century. 
The question that was addressed was posed by George Darwin: 

"The subject of the fixidity or mobility of the earth's axis 
of rotation in that body... [has] from time to time attracted 
the notice of mathematicians and geologists. The latter 
look anxiously for some grand cause capable of producing 
such an enormous effect as the glacial period. Impressed 
by the magnitude of the phenomenon, several geologists 
have postulated... a wide variability in the position of 
the poles on the earth; and this, again, they have sought 
to refer back to the upheaval and subsidence of continents."^ 

Very minor movement of the poles of perhaps 50 feet is measured 
by modern astronomers and is referred to as the variation of latitude. 
This variation is thought to be due to seasonal changes in the distribu 
tion of the air masses over the earth. 

^G.H. Darwin, Scientific Papers , (Cambridge, 1910), Volume 3, p. 1. 
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The use of floating frames was introduced in order to formulate 
equations of motion for a body undergoing small deformations due to 
elevation and subsidence of continents and sea beds. The researches 
of George Darwin were the most detailed in the evaluation of specific 
results. He introduced a principal axis frame and derived equations of 
motion for this frame. By adopting a system consisting of a rigid 
earth and a deformable thin shell, he studied the effects of deforma- 
tion. The shell Xi7as deformed in such a way as to give the maximum 
movement of the poles consistent with geological evidence. He was able 
to show that a movement of the poles of 8 degrees would require one 
half of the earths surface to be deflected by nearly two miles. This 
would effectively make continents out of oceans and vice versa. Darwin 
then concludes: 

”lf the geologists are right in supposing that where the 

continents now stand they have always stood, would it 

not be almost necessary to give up any hypothesis which 

2 

involved a very wide excursion of the poles?” 

The mean axes which were popularized by Tisserand and which bear 
his name were introduced by Gylden in his study of the rotation of the 
earth. This is stated by Tisserand in Sections 214 and 216 of his 
Traite de Mechanique Celeste. Tisserand derives the equations of 
motion for what we call the Tisserand frame by setting the internal 
angular momentum to zero. He reports then on the work of Glyden in 
using this frame in deriving an expression for the deviation of the 
mean axes from the principal axes of the system. The application of 

^Ibid . , p . 39. 


24 


the floating reference frame was similar in intent to that of Darwin, 
but his conclusions regarding the motion of the poles was more analyti- 
cal in nature and not as concrete. 
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CHAPTER 3 


THE TISSERAND FRAJIE FOR A DEFORMABLE BODY 
3. 1 Chapter Summary 

This chapter focuses in detail on the use of the Tisserand frame 
for vibration analysis of unconstrained systems. The choice of mode 
shapes to be used for expansion of the deformations relative to the 
floating frame is crucial. With an unwise choice of mode shapes, the 
constraint relationship may depend on all of the coordinates and may be 
difficult to deal with. The simplest constraint relationships result 
from the use of the free-free modes of an unconstrained system. The 
orthogonality of the rigid body modes (zero natural frequency) to the 
deformational modes allows for a rather elegant ^ priori evaluation of 
the constraint relationships. The constraint relationships then 
Involve only the rigid body modal coordinates, which must be zero. 

This result may then be interpreted as a requirement that the rigid 
body mode is fixed relative to the Tisserand frame. The remaining 
coordinates are then independent, a condition which can greatly simplify 
the formulation of equations of motion. This special case of the 
Tisserand constraint that uses the free-free modes of an unconstrained 
system is termed the Mode Shape Constraint. An example is worked out 
in order to demonstrate the properties of the Mode Shape Constraint. 

The classical Tisserand constraint cannot be employed for systems 
with rotating internal members. This class of systems represents space- 
craft with momentum exchange controllers. These represent reaction 
wheel, momentum wheel and control moment gyro control systems. Several 
alternate extensions to the Tisserand frame are discussed and a Modified 
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Tisserand Constraint is introduced. The properties of this constraint 
relationship are discussed and the applicability of the Mode Shape 
Constraint is demonstrated. 

3 . 2 Mode Shape Constraint 

There were two requirements made in order to define the Tisserand 
frame for a deformable body. The first was that the center of mass be 
fixed in the frame. This would then insure that the system possessed no 
linear momentum relative to the frame. But by introducing this con- 
straint* three additional variables were required to specify the loca- 
tion of the point of the frame occupied by the mass center, which point 
may be referred to as the frame origin. The second defining require- 
ment was that the system had no angular momentum relative to the frame. 
This constraint then introduced an additional three coordinates to 
specify the angular orientation of the frame. The net result of 
injecting the Tisserand frame into the problem for a deformable body 
is to expand the dimension of the problem by six and to interrelate 
the coordinates by six scalar equations (two vector equations). It is 
the interrelation of the coordinates, and the resulting fact that they 
are no longer linearly independent, which can prove troublesome in 
the formulation of the equations of motion. Constraint relationships 
require the use of Lagrange multipliers when Lagrange's equations are 
used. It would be extremely advantageous to evaluate the constraint 
relationships and reduce the order of the system before the equations 
of motion are formulated. 

There exists a very intimate relationship between the Tisserand 
frame and the free-free modes of an unconstrained system; this special 
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relationship we call the Mode Shape Constraint. This relationship 
employs the orthogonality conditions that exist between the rigid body 
modes and the deformational modes (which have nonzero natural frequen- 
cies) to allow a rather elegant ^ priori evaluation of the Tisserand 
constraint. This is accomplished by setting all rigid body modal coor- 
dinates to zero. This effectively reduces the order of the system to 
the original value and the remaining coordinates are all independent. 

The constraint relationships that define a Tisserand constraint 
for a deformable body (which by virtue of small deformations are also 
the Buckens constraints) are given by 


I 


p dm = 0 


( 3 . 1 ) 


and 


•'D 


dm = 0 


( 3 . 2 ) 


where £ and £ are now measured from the frame origin to the generic 
mass element. A separation of variables is introduced using the free- 
free mode shapes of an unconstrained system 
n 


u . . 

j=0 J 


( 3 . 3 ) 


where the mode shapes for distinct eigenvalues are orthogonal 

dm = 0 i 9 ^ j . (3.4) 

D 

The constraint relationships are then written as 


ti 

1 = 0 * 1 ) 


dm Oj = 


( 3 . 5 ) 
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and 


2 / £ dm n = 0 (3.6) 

j=0 •'d J 


where £ is the undeformed location of the generic mass element. In the 
above, the zero subscript has been used to refer to all of the rigid 
body modes. These consist of a total of six modes. The three transla- 
tional rigid body modes are taVien as uniform translation along each of 
the axes. They mav be written as 

(3.7) 

The three rotational rigid body modes are taken 
as small rotations about each of the axes. They are pictured in Figure 
3.1 and are written as 




0 

ii = 

0 

*2 ■ 

0 

ia ' "£a 


where C is a constant. 


,0 y a z £ 

■^4 fi , -3 “ £ -2 

.0 z a X 2 

■^5 £ -1 " £ -3 

,0 X a y a 

£ -2 £ -1 


(3.8) 


The translational rigid body modes are orthogonal to the rotational 
rigid body modes since the undeformed position has its origin at the 
center of mass and 




p dm = 0 . 



(3.9) 


It will be most fruitful in the evaluation of the constraint relation- 
ships to investigate the orthogonality properties of the free-free 
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inodes further. Let be a deformational mode corresponding to a 
nonzero natural frequency. 




= {f}" 




(3.10) 


This mode will be orthogonal to the translational modes, yielding the 
relationships 


/ ^9 * ^ dm = C / (|)^ dm = 0 

•'d ^ 

^ ^2 * -^ dm = 0 


X (J>9 • (p. dm - C I <})^ dm *= 0 

-3 '^3 

This may be \<rritten in the vector form as 


(3.11) 


/ (p^ dm = 0 


J 


(3.12) 


'D 


This is recognized as the coefficient in the center of mass constraint 
of n . Coupled with the orthogonality of the rigid body translational 
modes to the rotational modes shown in Equation (3.9), the center of 
mass constraint can be rewritten as a relationship Involving only the 
translational rigid body modal coordinates 

0 


^ 1 ° ^01 "^02 '’03 “ “ 
Evaluating the integrals will yield the three scalar equations 
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(3.14) 


Ifll ' 0 
•^102 " ® 
u(rno3 = 0 

where %M is the total system mass. It is now quite clear that the 
center of mass constraint requires that the translational rigid body 
modal coordinates be zero 

'loi ■ “ 

%2 - “ 


'^03 ® 


(3.15) 


Now we shall turn to the constraint setting the internal angular 
momentum to zero by studying the orthogonality of a deforma tional mode. 

to the rotational rigid body modes. The following orthogonality 
conditions hold 

dm = 0 




dm = 0 


dm = 0 


C3.16) 


When the expression for the rigid body modes (Equation (3.8)) 
substituted into the above, we have 


are 


z(|)^) dm = 0 


Ayri - 

: -^1. . 

jf(z(|)^l- x<l)^) dm = 0 
jf(x{{)^ - y(f>^) dm = 0 


(3.17) 
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vrhere the factor of Z inverse has been eliminated. If the three scalar 
equations are considered as components of the vector bases of the float- 
ing frame, then one vector cross product relationship results from 
Equation (3,16) 


/ £ X ({)- 


p X dm = 0 


where p is the position of the undeformed generic mass element 


(3.18) 


p - {frt 


(3.19) 


This vector cross product expression (Equation (3.18)) is recognized as 
the coefficient of the modal coordinate velocity rij in the constraint 
relationship (Equation (3,6)). The constraint relationship can then be 
written without the deformational modal velocities and also without the 
translational rigid body modal velocities. These latter are set to 
zero as a result of the center of mass constraint (Equation (3.15)), 
Writing the simplified expression for the constraint gives 

^ I 4 +J* £ >^ ^5 dm + J*£ X dm = 0 (3.20) 

Substituting the expression for the rigid body modes into the above 


yields 


{f}^ 


m 

i- 

iu- 


(yVz^)dm^ 


y X dm I riQ^ + 


z X dm 


li- 

\L- 


-X y dm 


(x^+z^)dm 


^05 


z y dm 


/ 


■ “ \ 
l- 
U 


y 2; dm 


2 2 
(x +y 


'06 


V= 0 


(3.21) 
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The moments and products of Inertia of the undeformed system are readily 
identified in the above. Writing the resultant scalar equations in 
matrix form yields 






^ - 

"ll 

’•12 

’'13 


^04 

^21 

.’'22 

^23 


^05 ■ 

J-31 

’■32 

’'33_ 


."^06. 


Since the inertia matrix is positive definite, 
this constraint equation is 


(3.22) 

the only solution of 





(3.23) 


This evaluation of the constraint relationship using free-free modes 
states definitively that the rigid body mode frame is fixed relative 
to the Tisserand frame. If initially aligned, the two frames will 
remain coincident. This is the prime advantage of the Mode Shape Con- 
straint and is the connection between the rigid body mode and the 
Tisserand frame that was mentioned at the end of Section 2.3. Thus the 
Mode Shape Constraint, which involves a Tisserand frame, can be used to 
define a rigid body mode frame. This result has been discussed by 
Likins,"^ but has been expanded upon here. The mention of the connection 
between the orthogonality of the rigid body and deformational modes and 
the momentum expression x/as made by Buckens. He used it only to set 


Peter W. Likins , "Analytical Dynamics and Nonrigid Spacecraft Simula- 
tion" (JPL TR 32-1593) p. 22. 
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total momentum to zero. Although the frame introduced by Buckens is 

crucial in the above arguments since it neglects second order terms, 

Buckens did not use the orthogpnaiity properties in order to evaluate 

2 

the constraint relationship for ai floating frame. 

3.3 Discussion and Example of the Mode Shape Constraint 

The development of the Mode Shape Constraint in the previous 
section has been quite detailed and may seem complex. A discussion 
and summary of the pertinent results are in order. The example problem 
from Section 2.4 will be expanded in order to demonstrate the Mode 
Shape Constraint in a concrete manner. 

The Mode Shape Constraint is defined as a specific case of the 
Tisserand constraint where the deformations are expanded in terms of 
free-free modes of the unconstrained system. The main property of the 
Mode Shape Constraint is that it eliminates the rigid body modal coor- 
dinates from the problem and leaves the remaining coordinates indepen- 
dent. If an arbitrary choice of mode shapes is made, then the Tisserand 
constraint will render the modal coordinates dependent on each other. 

The specific structure of the interrelation of the modal coordinates 
depends on the choice of mode shapes. This may be determined from 
examining Equation (2.21). 

The procedure that results from applying the Mode Shape Constraint 

is to expand the deformations in terms of the deformational modes only, 

ignoring the rigid body modes, which have zero natural frequency. The 

constraint relationship is fulfilled because each coefficient in 
_ 

F. Buckens "The Influence of Elastic Components on the Attitude Stab- 
ility of a Satellite," Proceedings of the Fifth International 
Symposium on Space Technology and Science (Tokyo, 1964), p. 196. 
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Equation (2,21) is set to zero. The coordinates are then Independent 
because of the specific choice of mode shapes. This interpretation 
gives rise to the name "Mode Shape" Constraint. The Mode Shape Con- 
straint may also be viewed as a straightforward way of working with 
the rigid body mode frame. As noted before, the constraint relation- 
ship requires that the rigid body modes do not move relative to the 
Tisserand frame. Thus, the mode shape constraint may be used to locate 
both frames. This is a much easier method to locate the rigid body 
mode frame than a requirement involving zero strain energy. 

The discussion will now take up the example problem first intro- 
duced in Section 2.4. This is a case of planar motion of a long slender 
member, so the constraint relationship will simplify to the scalar form 


I 


X u dm = 0 . 


(3.24) 


If the velocity is expanded in terms of mode shapes, then 


t(J[. = 


0 


(3.25) 


The mode shapes for the system shown in Figures 2.1 and 2,2 are the 
rigid body mode and one deformational mode relative to the Tisserand 
frame. They are shown in Figure 3.2. The mode shapes are 


iO X 

t’ “I 


- .586 X 
— r + .414 X 


A. 




- <x < 0 

V2 

0 < x < £ 


(3.26) 
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Several characteristics of the mode shapes must be mentioned. 

They represent only small deformations and always place the origin at 
the center of mass (Equation (2.24)). The deformational mode resulted 
from setting the relative deformation, a, to one. The rotation of the 
body fixed frame relative to the Tisserand frame, 3» would then be .586 
(Equation (2.41)). The resultant movement of the center of mass would 
be - ^ . The essential property of these two modes is that they are 
orthogonal. 





(3.27) 


Thus the mass matrix is diagonal. Since the relative deformations of 
the system result from only the second mode shape, the stiffness matrix 
is zero except for ^ 22 ' modes chosen therefore are the natural 

modes of the unconstrained system which we call the free-free modes. 

The system is unconstrained as the singular stiffness matrix would 
indicate. 

Now let us use the free-free modes in order to evaluate the 
coefficients of the Tisserand constraint relationship shown in Equation 
(3.25). 





= 2.218 m 

X (j)^ dm = 0 (3.28) 

The last coefficient is set to zero because of the orthogonality condi- 
tion. The Tisserand constraint for this system may be written 
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zero. 


Hq + X dmj = 0 

(2.218 m £^)nQ= 0 (3.2^) 

The constraint requires that the rigid body modal velocity be 
Thus the rigid body mode does not move relative to the Tisserand 
frame. Since the original choice of axes is arbitrary, the modal coor- 
dinate may be set to zero and will remain zero because of the constraint. 
3.4 Modified Tisserand Constraint 

This section will deal with an extension of the Tisserand frame 
for deformable bodies to cover a new class of systems. Consider a 
deformable body with small relative displacements to which is added a 
spinning rigid rotor. The new system will no longer be a strictly 
defoirmable body, since the movement of the rotor will Involve large 
relative displacements. If the classical Tisserand constraint is 
applied to this system, the frame itself must rotate relative to the 
system. The system angular momentum will consist of a component from 
the deformable body and a component from the rotor. The frame must 
rotate in such a way as to carry the total system angular momentum and 
thereby set internal angular momentum to zero. The frame motion cannot 
coincide with deformable body component of the system and large deforma- 
tions result. 

Two methods of modifying the Tisserand constraint present them- 
selves. The first is to define a Tisserand frame for the deformable 
body only. This approach ignores the rotor altogether when the internal 
angular momentum relative to the deformable body center of mass is set 
to zero. The constraint relationship is then written 
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(3.30) 



X u dm = 0 


when p_ is the undeformed location relative to the center of mass of 
the deformable body. This approach is not used because it does not 
achieve the greatest simplification. 

The second method to be considered is to examine the expression 
for the system angular momentum. The constraint relationship can then 
be made to simplify the structure of the angular momentum as much as 
possible while still assuring small deformations. 

Before we write the expression for the system angular momentum, 
let us define several vector bases: 

{_i}: Inertially fixed reference axes 

/\ 

{jf}: Axes fixed in the floating frame 

/N 

{_b}: Axer. fixed in the mass element at the rotor mounting point 

{v}: Axes fixed in the rotor 

The angular momentum of the system about its center of mass is 

H = / £ X £ dm (3.31) 

where S represents integration over the entire system, both the rotor 
and the deformable body. Evaluating the angular momentum for the 
deformable body will give 



X p dm 




/* p X u 

‘'d “ ~ 


dm 


(3.32) 


where is the deformable body dyadic about the system center of mass. 
For the rotor, static balance has been assumed. The similar expression 
integrated over the rotor, R, is 
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(3.33) 


_^£X^dm= []= • «“ + Dr • + «r Pg 1 r 


where 

Rotor inertia dyadic about the system center of mass 
Rotor inertia dyadic about its own center of mass 
Rotor mass 

Pj^.v Location of rotor center of mass relative to system center 
of mass. 

By combining these two quantities, the system angular momentum is found 
to be 


H 


= n * Dr * y*£ ^ \ pR ^ £r 


(3.34) 


The total system inertia dyadic is given by 


□ - Od + D 


(3.35) 


The total angular momentum expression can be significantly simplified. 
In order to do this, we introduce a Modified Tisserand Constraint. 


^£ ^ if. 3m + X = 0 


(3.36) 


Here we have made use of the chain rule for angular velocities 

+ u'’* . (3.37) 

The system angular momentum then becomes 

I = D * . (3.38) 

The above has the same structure as a rigid body with an attached rotor 
The essential difference is that the Inertia dyadic is not constant and 
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that the relative angular momentum of the rotor (the second term) is 
affected by the deformations. It is important to understand that rota- 
tions of the axes at the mounting point will reorient the rotor axes. 

The internal angular momentum of the system relative to the frame 
is no longer zero. 


p X u dm = . (3.39) 

This is a very simple result and may be verified by direct computation. 
Alternatively, a quick comparison of the above angular momentum expres- 
sion (Equation (3.38)) and the general angular momentum expression 
(Equation (2.3)) will give the same result. The one term evaluation of 
the internal angular momentum will facilitate the evaluation of the 
equations of motion. In order to demonstrate this, let us evaluate the 
kinetic energy. We start with the general result given above (Equation 
(2.6)) vzhich is now evaluated over the entire system 



T = f • £ dm + • [] • • y*£ X £ dm (3.40) 

The first term must be evaluated over the deformable body and the rotor. 

l/*o o, l/*o o 
-2 ^ £ • £ dm = -2 y u • u dm 


1 To o 

2 


1 ..vb 


2 “r £-E • -^E 2 • 


De ■ 


w 


vb 


+ • Dr • ‘ D„ • (3.41) 


- 2 - — R 

Making the substitutions indicated will yield the kinetic energy for a 
deformable body and attached rigid rotor with static balance when the 
Modified Tisserand Constraint is employed. 
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bf 


T = u • u dm + y • [] • 

D 

. 1 ,, o o 1 vb f~| vb . vb f~] bf , 1 bf n 

+ 2 ”r£r • -aR"^ lii ■ Ur • ii +“ • Ur • w +2“ 'UR-ii! 

(3.42) 

Even with the simplification that comes from the Modified Tisserand 
Constraint, this system has a complicated structure for the kinetic 
energy. 

The Modified Tisserand frame is distinct from the frame created 
by applying the classical Tisserand constraint to the deformable body 
only (Equation (3.30)). In order to demonstrate this, the internal 
angular momentum of the deformable body relative to its own center of 
mass will be calculated for the Modified Tisserand Constraint. Starting 
with Equation (3.39) we may write 


p_ dm p_ X dm = 


(3.43) 


Each component of internal angular momentum may be evaluated by rewrit- 
ing the expression relative to its own center of mass rather than the 
system center of mass. 


^ ^ X u dm + - Mj,)r^ ^ Lq 


f 

I p X p dm 
- 


=D.- 


vb , 

R U + 


Dr • + Mr X £r 


(3.44) 


Where Is the system mass and is the location of deformable body 
center of mass relative to the system center of mass. By using the 
center of mass expression 

(^- + Mg Pp_ = 0 (3.45) 
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we may evaluate the internal angular momentum of the deformable body 
relative to its oxm center of mass 


L 




X u dm 


= -D 


R - 


bf 



(3.46) 


This quantity will not generally be zero, and the deformable body will 
have a net internal angular momentum relative to its own center of mass 
when the Modified Tisserand Constraint is employed. This result also 
shows that the Modified Tisserand Constraint involves an extension of 
the classical Tisserand constraint for deformable bodies. 

Going back to total system angular momentum (Equation (3.38)) it 
will be Interesting to examine the frame motion. For the case where 
the system has zero net angular momentum and is free of external 
moments, the relationship will hold whereby 

D ‘ + Dr • 0 • (3.47) 

This relationship does not require that the frame be inertially fixed. 
Instead, the frame will move when the internal angular momentum (the 
second term) is altered in value. This is the process of momentum 
exchange that is employed by control systems. 

The most Important aspect of the Modified Tisserand Constraint 
has been left till now. We shall now attempt to tie together the 
advantages of the Mode Shape Constraint with the extended systems 
covered by the Modified Tisserand Constraint. The crux of this endeavor 
lies with the physical interpretation of the Modified Tisserand Con- 
straint (Equation (3 .36) ). Let us start by discussing the terms present 
in the expression. The first term is the internal angular momentum of 
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the deformable body. It Is the last two terms which Involve the rotor. 
But neither term involves the rotor spin rate. Rather, these two terms 
represent the Internal angular momentum of the rotor if it were indeed 
not spinning relative to the axes fixed at the mounting point. That is, 
the rotor would be fixed in the reference frame b. The constraint 
relationship may now be interpreted as a classical Tisserand constraint 
for a system with a "frozen" rotor. Once this Interpretation has been 
made, the Mode Shape Constraint associated with a classical Tisserand 
frame may be applied to a deformable body with an attached rotor. The 
mode shapes used to expand the deformations are found from the free- 
free analysis of the system where the rotors have been frozen. Thus 
the rotor mass and inertia properties contribute to the free-free modes. 
The Modified Tisserand constraint will now require only that the rigid- 
body modal velocities be zero. This allows them to be ignored, and the 
analysis proceeds with the frame variables and the modal coordinates 
associated with nonzero natural frequencies. All of these coordinates 
are independent since the constraint relationship was evaluated by 
suppressing the rigid body modes. 

All of the relationships given above for a system with one rotor 
with static balance will hold if additional rotors are added. The only 
change is that expressions involving rotor quantities are summed for 
all rotors . 
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CHAPTER 4 


FINITE ELEMENT EQUATIONS OF MOTION 
4.1 Chapter Summary 

Thl&v chapter deals with the derivation of a generalized set of 
equations of r>otion for the class of spacecraft with distributed 
flexibility and momentum exchange controllers. The Finite Element 
Model adopted for this discussion consists of a finite number of rigid 
bodies interconnected by massless, elastic elements. This model allows 
the equations of motion of a continuum to be represented as a finite 
number of ordinary differential equations instead of a partial differ- 
ential equation. The Finite Element approach is particularly powerful 
for dealing with systems characterized by nonuniform mass and stiffness 
properties. The Finite Element Model introduces six coordinates for 
each rigid body in the model, three coordinates for translation and 
three for rotation. The total number of discrete coordinates is then 
seen to be large even for a system of modest complexity. This suggests 
the adoption of distributed coordinates obtained by a modal analysis 
and the truncation to a smaller number of modes . The process of trunca- 
tion is complicated and will be discussed in Chapter 6. 

The equations of motion developed here deal with a system which 
has small, first order displacements and rotations, relative to an 
Inertlally fixed coordinate system. This requires that the variation 
of the rotor angular momentum is also a small, first order quantity. 
Although the equations of motion for this system can be written using 
an Inertially fixed coordinate system, a Modified Tisserand Frame will 
be introduced. This derivation will demonstrate the process of 
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formulating equations of motion when this floating frame is employed, 
as may be necessary for the case of large system rotations. Particular 
attention will be paid to the impact of the constraint relation. 

The equations of motion will be formulated using Lagraxige’s equa- 
tions. The first step consists of the calculation of the system 
kinetic energy. In order to arrive at equations of motion containing 
first order quantities, care must be taken to include all second order 
quantities in the kinetic energy. Because the rotor angular velocity 
is a large quantity, not all terms can be represented by their first 
order approximation. The generalized forces for the system are conser- 
vative and the strain energy is a simple quadratic form involving the 
stiffness matrix. 

The use of the Modified Tisserand Frame will involve a constraint 
relationship. Because the generalized coordinates are not independent, 
Lagrange multipliers are introduced in the equations of motion. The 
modal analysis is performed for the free-free system with the rotors 
frozen and the transformation from discrete to distributed coordinates 
does not include the rigid body modes. Thus, the Mode Shape Constraint 
is Introduced td the problem and the result is the elimination of the 
Lagrange multipliers. The final set of equations of motion involves 
coupling between vibration and the Tisserand Frame motion. Gyroscopic 
coupling between the vibration coordinates is also present. If damping 
is desired in the system, a term may be added in the final set of equa- 
tions to represent modal damping. 

This set of equations could then be truncated to the frame 
variables and a smaller number of necessary modal coordinates. These 
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equations would be representative of the system motion of a spacecraft 
with small displacements from inertially fixed axes. They could then 
be incorporated in a control system design to attain precise pointing. 

4.2 System Kinetic Energy 

The first step in the process of formulating the equations of 
motion using Lagrange’s equations is the derivation of the expression 
representing the system kinetic energy. Several assumptions have been 
made to facilitate this endeavor. The system is assumed to be free of 
external forces and the system center of mass then represents an iner- 
tial reference point. The resulting equations of motion will therefore 
not include any expression representing accelerations of the system 
center or mass . The equations of motion representing translation of 
the system center of mass will not appear in the derivation since in 
effect they have been solved ^ priori . The assumption is made that 
the rotor center of mass is located at the center of mass of the 
attached body. This assumption simplifies the structure of the result- 
ing equations, but could be abandoned easily if it were judged inade- 
quate for analysis of the system at hand. The rotor is assumed to 
be statically balanced. The angular momentum of the rotor is also 
assumed to undergo only small, first order variations. This would 
correspond to a requirement that the rotor produces only small, first 
order torques. This Implies only first order changes in rotor speed 
for a momentum wheel or only first order changes in gimbal angles for 
a control moment gyro (CMG) . A Modified Tisserand Frame is introduced 
for the system. This results in a significant simplification of the 
internal angular momentum expression Xi?hich appears in the kinetic energy. 


49 



As shown before (Section 3.4), the expression for the kinetic 
energy for a deformable body containing a spinning rotor is 


T = 


1 I o o f o 


dm 




D • “ 


fi 


(4.1) 


If a Modified Tisserand Constraint is introduced to define the floating 
frame, then the constraint is written as follows for a rotor with 
static balance. 


^£X£dm + Mj^Pj^X£^ + []^ • ^ = 


0 


(4.2) 


The resulting expression for the Internal angular momentum is 

J X ^ dm - • 0^^ 

The kinetic energy can now be put in the form 

R 


(4.3) 


T='|y*£‘£dm + 




D • • 


(4.4) 


In order to analyze the first expression in the kinetic energy, the 

location of the generic mass element relative to the system mass center, 

£, is defined for each body as in Figure 4.1 

^ -i . i . i 
p e p + u + r 


(4.5) 

where £^ is the undeformed position of center of mass of body i and is 
fixed in the frame f . 


—1 

p = 0 . 


(4.6) 
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The displacement In frame f of the center of mass of body 1 is 
represented by the vector, }^ . The location of a generic mass element 
of body i relative to the body center of mass is ^ and therefore obeys 
the constraint 



r dm = 0 


(4.7) 


where the integration is taken over body i. A similar expression 
results for the location of a generic mass element of the rotor 


-R ^ 

p = p + 


R R 
u + r 


(4.8) 


where the superscript R refers to the rotor. The derivatives of the 
mass element expressions relative to the floating frame give 


°1 oi b f _ i 
p = u + w X r 


°R OR . vf ^ R 
p = u + CO X r 


(4.9) 


where the reference frame b is fixed in the body i, the reference frame 
V is fixed in the rotor and f is fixed in the modified Tisserand frame. 
Remembering the center of mass constraint (Equation (4.7)) and the chain 
rule for angular velocities 






(4.10) 


we may now write the vector-dyadic expression for the system kinetic 
energy 
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^/i oi oi 1 b^f 

= 2 - ■ U H / 

1 °R .°R , 1 ,.vb^ . 

T ^ "2 — LI ’ — 


^ 1 , b^f riR , ^ 

■'■ "2 — * LI * — — 




vb 


R 


+ 0) 


fi 


D 


R 


, vb^ . 1, fi 
— "2 — 


D-“ 


fi 


(4.11) 


R . 


In the above expression, the reference frarae b is fixed in the body to 

which the rotor is attached and the inertia dyadics and are for 

body i and the rotor, respectively, about their individual centers of 

mass. The body i and rotor masses are m and M respectively. The 

i K. 

number of bodies in the model is n, not counting the rotor. 

The next step in the evaluation of the system kinetic energy is 
to adopt a set of vector bases and express the kinetic energy in matrix 
form. The following set of vector bases are adopted: 

A 

{O; Fixed orientation relative to the inertially fixed reference 
frame. 

A 

{f}: Fixed orientation relative to the Modified Tisserand frame. 

{_b^}: Fixed orientation relative to the locally attached frame 
in body i. Superscript R refers to the body to which the 
rotor is attached. 

/N . ■ ■ 

{jv}: Fixed orientation relative to the rotor attached frame. 

Care must be taken in representing the vectors present in Equation 

(4.11) in terms of first order approximations. This cannot be done for 

vb^ 


the following terms since to 


0) 




is not a first order quantity 
R 


□ R vb 

• (£ 


(4.12) 
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b^f f i 

The angular velocities, co and Oi , present here must include second 

order terms in order to keep all the second order terms present in the 
kinetic energy. Referring to Appendix B, these terms may be written 

0^ ^ = {bV({3^} + [3^] {s’"}) 


0) 


0) 


= {f}^(fe} + [0^] {0}) 


(4.13) 

The other vectors in the kinetic energy may be represented by their 
first order forms 

u" ={f}^{u"} 

= {b^}^ 

= {b^}^ [I^] (bH 
= {f}"^ {u^} 




[]^ - rCRlT , rCRi 


= (ry [i^] & 


[] = {f}^ [I] {?} 


(4.14) 

The matrix expression for the system kinetic energy may now be written 
in a form which includes all second order terms 




+ -i + i + i {s}’^[i^]{s} 


R^ 


+ ]{s} + + {8}’^[P][T Hs} 


R 


R^ 


R, 


+ {6}'^[ijji{s} + {ef [ej[i„]{s) + -i {fi}''m{e} 




(4.15) 
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In the above expression the "tilde" operator represents a skew symmetric 
matrix as typified by 


[3^] - 


0 


-^3 


0 - 3 ^ 

0 


,R 


(4.16) 


The requirement that the variation of the rotor angular momentum be a 
small, first order quantity may now be applied to the components repre- 
sented in Equation (4.14). The rotor angular momentum is written 


(4.17) 


Since the angular velocity relative to inertial space includes the 
first order terms m and a they can only contribute first order 
quantities to the rotor angular momentum. The variation of the rotor 
angular momentum which must be kept first order is: 

, R 


If ('V <=>) 


(4.18) 


For a reaction wheel or momentum wheel, the inertia matrix will be 
constant and the first order requirement will be 

(4.19) 

For a control moment gyro (CMG) , the rotor speed will be constant. 

Then the first order requirement will be 

T ■ ■ 


'V ■ 0^ 


(4.20) 


These first order requirements will prove Important when the equations 

of motion are derived using the kinetic energy. 
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4.3 Constraint Equations 


Since the kinetic energy was evaluated In terms of a Modified 
Tlsserand Frame, the next step In the formulation of equations of 
motion Is the evaluation of the ssaoclated constraint relationship. 
As given before (Equation (4.2)), the Modified Tlsserand Constraint 
for a deformable body and a rotor with static balance Is 


£ X 5 dm + X 5^ + - 0 


(4.21) 


If the same conventions are adopted as In the previous section (Figure 
5.1), the first term may be evaluated by using Equations (4.5) and 
(4.9). The resulting first order expression is then 

^ ^ Dr ' ° 

The vectors and dyadics may be expressed in the same manner as 
before (Equations (4.13) and (4.14)) when a system of vector bases Is 
adopted. To this representation we must add 
* {f}’^ 

= {f}’^ {x'^} (4.23) 

A common vector basis may be adopted, because to the first order 

{b^}^ = {i}^([E] + [3^]) (4.24) 

where [E] Is an identity matrix and the "tilde" is the skew-s 3 Tmnetric 
representation shown in Equation (4.16). The resulting three scalar 
constraint equations can be written in the form 

+ [l'^]{3^>)+ M^[#](u^} + [Ij^]{3^} = 0 (4.25) 

Here the "tilde" skew-symmetric operator is a representation of the 
cross product . 
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4.4 Equations of Motion-Discrete Coordinates 


The equations of motion may now be formulated using Lagrange's 
equations. For this discussion, the generalized coordinates are not 
independent and Lagrange multipliers must be introduced. The form of 
Lagrange’s equations to be used here is 


A_ 

dt 




Sk 

k = 1, ... , V 


(4.26) 


where V is the number of generalized coordinates. There are m con- 
straint relationships of the form 


^Sk ° 

k=l 

S = 1, . . . , m (4.27) 

For our problem, V is equal to 6n + 3, and m is equal to 3. The resul- 
tant equations of motion are first order only. It is important to 
realize that the application of Lagrange's equations to the system 
kinetic energy (Equation (4.15)) yields some second order terms. A 
crucial part of the equation derivation process is the identification 
and elimination of all second order terms; the requirement that the 
variation of the rotor angular momentum be first order will identify 
many second order terms to be eliminated by recognizing as first order 
quantities the terms in Equations (4.19) and (4.20). 

The equations will be derived in groups of three. The first set 

of equations is for the Tisserand frame variables, {0}. As shown in 

fi 

Appendix B, the second order term in the angular velocity, W , can 
be written in two forms to facilitate the taking of derivatives 
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(4.28) 


{6)’'[0^)'' = {e}’^ [63]'' 

The second order terms have a very interesting structu? ii“ which will be 
used here 

- tejl’' = I«1 . (4.29) 

Once the above are recognized and all second order terms are eliminated, 
the system rotation equations (for rotation of the Modified Tisserand 
Frame) are written 

[iHe) - ([i^Hs>){e} - 

+ [Ir]{s> + [ipj{s} = {T^} (4.30) 

where {Tg} represents externally applied torques. 

The next sets of equations will deal with the translation and 
rotation of body i. This body represents all bodies in the model which 
do not have attached rotors. The translation equations are 

m^{u^} = + m^[x^]0(} . (4.31) 

Here the coordinates are present in the constraint relationship and the 
three Lagrange multipliers are introduced in a column matrix {A}. The 
generalized forces are represented in a column matrix {Q^}. The rota- 
tion equations for body i are 

= {Q^} - (4.32) 

Equations will now be written for body R, which does have an attached 
rotor. The translation equations are 

(itij^ + = {Q^} + (nij^ + Mj^)[S^]{X} (4.33) 
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The rotation equations for body R are more complicated. Use must be 
made of equations similar to Equations (4.28) and (4.29). When all 
second order terms have been eliminated, these equations can be written 


([i"] + [IrIH?) - ([ipHs}){Sh + [ig]{s} + [y{s) 


- (ti^Hs}){e} - {q„} - ([i"]’' + [I„]ba) 


R-' 


(4.34) 


The inertia matrix of body R about its own center of mass is [I^] and 
the rotor inertia matrix about the same point is 

R 

The derivation of the equations of motion is now complete. In 
order to examine their Structure, a column matrix {q} will be introduced, 
This represents the translation and rotation displacements for all 
bodies . 

/{u^} 


{q} = 




Kb^} 


(4.35) 


{u"} 

{B^} 

The generalized forces for all the bodies may be rewritten in terms of 
a stiffness matrix, [K] . 


{Q> = 




I{q"} 


= - [K]{q) 


(4.36) 
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It is now possible to put the equations describing the translation and 
rotations of the bodies in a compact form which we will call the system 
vibration equations 

[M]{q> + [G]{q> + [K]{q} = [R]{e> + {F} + [L]{X> (4.37) 

The component matrices in the system vibration equations have the form ' h 




0 


[R] = 


[L] = 


I 

I 

0 

[IrHs} 

0 

I 

I 

I 

0 


(4.41) 


R, 


(t\ + M^) [x ] 


xujx"] 

- [l"]^ 


(4.42) 


The system iTotation equation may be rewritten using equations 
(4.35) and (4.41) . 


tiHS} - ([l^]{s}){0} - - tRl’^lq} - tlgUs} - [ijHs} (4.43) 

Here the external torques have been set to zero. The term involving 
the derivative of the rotor inertia matrix in the above equation and 
in Equation (4.40) can be expressed in an alternate way. 

= [s][Ij^]{s> (4.44) 

The complete equations of motion for the system have now been 
found. The next step is to attempt to introduce distributed coordinates 
into the analysis. 


6l 



4.5 Modal Analysis 


The system equations of motion have now been formulated using 
discrete coordinates and consist of the system rotation Equations (4.43), 
the system vibration Equations (4.37) and the constraint Equations 
(4.25). The modal analysis will be introduced in order to substitute 
distributed coordinates for discrete coordinates in the vibration equa- 
tions. The objective is to reduce the order of the system (number of 
degrees of freedom) while still being able to represent system behavior 
adequately. This process is generally termed "truncation*'. 

The first decision is the choice of the eigenvalue problem to be 
solved. If the eigenvalue analysis is performed on the homogeneous 
vibration equation 

[M]{q} + [G]{q> + [K]{q} = 0 (4.45) 

the result will yield complex modes, or eigenvectors represented by 
complex numbers. A major drawback of this course of action is that 
these modes depend on the spin rate. From a numerical standpoint, 
working x^ith complex numbers may also prove cumbersome. 

Real modes would result from an eigenvalue analysis of the equa- 
tions 

[M]{q> + [K]{q} = 0 . (4.46) 

since [M] and [K] are real and symmetric. This represents a modal 
analysis of the free-free system where the rotor has been "frozen" 

(spin rate set to zero). Besides the numerical ease of working with 
real numbers only, the major advantage of this choice of eigenvalue 
problem is that the Modified Tlsserand Frame may be defined by the Mode 
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Shape Constraint. This will eliminate the Lagrange multipliers from 
the vibration equations. 

Let [(()] be a rectangular modal matrix corresponding to Equation 
(4.46) where the six rigid body modes have been eliminated. The defor- 
mational modes included in [cf>] have nonzero natural frequencies and are 
normalized so that 


= [E] 

] 

where [E] is the identity matrix and 
natural frequencies squared. The substitution is made that 


[> 2 ] . 

0 ) 1 

. 


is a diagonal matrix of 


[Cj>] [K][(})] = 03 , 


= 03 


(4.47) 


{q} = (4.48) 

where {ri) depends on time only and the vibration equations are pre- 
multiplied by [(f)] transpose. The resulting vibration equations are 

{r\} + 

+ + [t|3]^{F} + [(|3]^[L]a} (4.49) 

An examination of the constraint relationship Equation (4.25), shows 
that it can be written in the matrix form 

- [L]^{q} = 0 (4.50) 

With the modal coordinates defined in Equation (4.48), the constraint 
relationship becomes 

- [L]^[(|)]{n} = 0 (4.51) 

But with the Mode Shape Constraint, the coordinates {f)} are independent 
as a consequence of the orthogonality of the nonzero frequency modes to 
the rigid body modes. Thus the constraint relationship must be 
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fulfilled by setting all coefficients of the modal coordinates to zero 
in the constraint equation 

- [hf[<^] = 0 (4.52) 

It is an easy step from the above to eliminate the Lagrahge multipliers 
from the problem by realizing 

[<})]^[L] = 0 . (4.53) 

When the modal matrix is introduced into the system rotation equation, 
the resulting equations for system behavior become 

[i]{e} - (ii^} {0} = - [A]'^{ri} - [i„]{s> - [ipHs} 


R- 


{?i} + 2 




(4.. 


{n} = - [<})]^[G]tc{)]{n} + [A]{0> + [cj)]'^{F} 


where the coupling between system rotation and vibration depends on 


[A] = HV[R] 


(4.55) 
f 

Modal damping is also added for this set of equations. The quantity 


lS 


is a diagonal matrix containing the percent of critical damping for 
each mode. This final set of equations describes the motion of a system 
which is free of external forces and torques and for which changes in 
the rotor angular momentum are first order. 

The structure of the system equations of motion has several 
Important aspects. First, the rotation of the Modified Tisserand Frame 
(rigid body mode) is coupled to the vibration and vice-versa. The 
second aspect is that the vibration modal coordinates are coupled to 
each other by the gyroscopic term. This complicates the truncation 
process. An examination of this gyroscopic term will have to be made 
in order to successfully truncate the number of modal coordinates. The 
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last point is that the "inertial coupling" present in the hybrid 
coordinate technique for analyzing flexible appendages is not present 
here. That is, the coupling is not between the second derivatives of 
the frame coordinates and modal coordinates. The use of the Modified 
Tisserand Frame has permitted the elimination of this "inertial 
coupling" from the final equations. 
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CHAPTER 5 


CONTINUUM ANALYSIS 

5 . 1 Chapter Suinmary 

In this chapter we define a very simple physical system for an 
analysis as a continuum. The continuum analysis uses Hamilton’s 
Principle to derive governing partial differential equations with 
boundary conditions. The physical system consists of a long slender 
beam with an axisymmetric rotor mounted at the tip of the beam with the 
axis of symmetry normal to the beam and parallel to a principal axis of 
the beam cross section, when the beam is undeformed. The rotor is sym- 
metric with static and dynamic balance. The center of mass of the rotor 
lies on the beam centroid in order to simplify the presentation. The 
rotor spins with a constant rate plus a small, first order speed varia- 
tion. The beam is uniform and inextensible and is modeled as a 
Bernoulli-Euler beam. The beam is allowed to vibrate in the plane 
defined by the beam and the rotor axis, in the orthogonal plane contain- 
ing the beam, and about the centroidal axis (torsional vibration). No 
damping is assumed and the system is free of external moments and forces. 

The governing partial differential equations have a very simple 
structure. Three separate equations result for the system. They 
represent In-p lane bending, out-of-plane bending and twist. All of 
the equations are homogeneous and are the familiar forms associated 
with the bending and torsion of a uniform bar. On the other hand, 
the boundary conditions for the equations have a very complicated 
structure. Because of the gyroscopic coupling provided by the rotor, 
the in-plane bending and torsion are coupled through the boundary 
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conditions. In addition to this coupling, the boundary conditions 
depend on the eigenvalues of the system. The out-of^-plane bending is 
not coupled to the other system motion, but the boundary conditions 
depend on the rotor speed changes and are therefore time dependent. 

The solution of an equation with time dependent boundary condi- 
tions involves a special transformation which for some cases renders 
the boundary conditions time independent. For a modified system, this 
solution technique is presented for out-of-plane bending. If a tradi- 
tional separation of variables is employed for the coupled bending and 
torsion, the independence of each coordinate will cause only one solu- 
tion to be admitted; this is the trivial solution corresponding to zero 
displacements. This result forces the separation of variables to adopt 
a new structure. The spatially dependent coordinates are distinct for 
bending and torsion, but they both have the same natural frequency. A 
difference in phase is allowed between bending and torsion. The result 
of this solution form involves complex mode shapes. The interpretation 
of these modes leads to the physical solution where in-plane bending and 
torsion are ninety degrees out of phase. The rotor tip is then seen to 
trace out an ellipse where one axis represents pure bending and the 
other axis represents pure torsion. 

The last section of the chapter provides numerical solutions for 
a particular choice of parameters. The solutions for different rotor 
speeds are presented to show the dependence of complex modes on the 
spin rate. The actual solution for the natural frequencies and mode 
shapes for specific systems allows greater physical insight into the 
dynamics of the system. 
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5.2 Governing Partial Differential Equations 


In this section, Hamilton's Principle will be used to derive the 
equations of motion and boundary conditions which govern the behavior 
of a very simple physical system. This system consists of a long 
slender flexible beam with a rotor mounted at the tip. This system is 
shown in Figure 5.1. The system is free of externally applied forces 
and moments. Because of this, the system center of mass is an inertial 
reference point. A set of inertially fixed unit vectors, {_i}, have their 
origin at the center of mass. The system has no damping and the beam is 
uniform and inextensible . It is modeled as a Bernoulli-Euler beam with 
degrees of freedom for displacement in-plane, out-of-plane and torsional 
twist about the centroidal axis. The beam does not include the effects 
of rotary inertia and shear deformation and no degree of freedom is 
included to allow longitudinal vibration. The rotor is mounted on the 
tip of the beam With the center of mass of the rotor placed on the 
beam centroid. The rotor is symmetric and has static and dynamic 
balance. A locally attached vector basis, {_b}, is fixed at the tip of 
the beam, with one axis aligned along the rotor axis of symmetry and 
another axis aligned along the beam. The origin of these axes is placed 
on the beam centroid. This set of axes is representative of a general 
set of axes fixed in each beam cross-section. A frame v is fixed 
relative to the rotor. The rotor speed consists of a constant rate 
plus a small, first order variation that is time dependent. It is not 
necessary to employ a Tisserand frame, since the requirement that the 
speed variation be small will ensure that the Inertial displacements of 
the beam are small (first order quantities). 
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A first order analysis is used when \/e apply Hamilton's Principle. 
For conservative, homonomic systems, Hamilton's Principle is written 



(5.1) 


where [£ is the Lagrangian of the system. The Lagrangian is the differ- 
ence between the system kinetic energy and the system potential energy. 
.2’ = T - U, (5.2) 


Hamilton's Principle requires that the actual path (Newtonian path) of 
the system in configuration space renders the value of the definite 
integral of Equation (5.1) stationary with respect to all arbitrary 
variations of the path between the two time points. The arbitrary 
variations from the true path create what is called the varied path. 

It must be understood in Hamilton's Principle that the varied path 
coincides with the true path at the endpoints of the integration; that 
is at time points t^ and t 2 the variations are zero. The advantage 
of Hamilton's Principle is that it provides in a formal manner both the 
partial differential equations of the system and the associated boundary 
conditions. 

In order to keep first order terms in the equations of .motion 
and the boundary conditions, second order terms must be preserved in 
the kinetic and strain energies. The kinetic energy for the system 
is first written as 



For the beam component of our system, the vector _p from the system mass 
center to the generic mass element may be written 
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£=r + £ = £+ u + £ (5.4) 

The geometry of the beam is shovm in Figure 5.2. The beam centroid is 
defined by the relationship for each cross-section dA 


£ 


£ dm = 0 


(5.5) 


The vector £ is normal to the centroidal axis and locates the mass 
element relative to the beam centroid. This vector is fixed in the 
locally attached vector basis {Jb} and has no component in the direc- 
tion. The vector is written 


£ = {b}^ 


(5.6) 


The torsional twist, is assumed to be about the centroid. If rotary 
inertia effects are to be ignored, then the only component of angular 
velocity of the locally attached frame, that needs to be considered 
xs that associated with torsion, ipb^. Thus the derivative of the vector 
£ with respect to inertial space is 

hj (5.7) 


bi 


£ = w ^ 2. s {^} 


The displacement of the beam centroid is represented by the vector 
jc. This vector may be broken into two components: a vector, £, which 

locates the position of the undeformed centroid and is fixed relative 
to the inertial axes, and a vector u. which represents the displacement 
of the centroid. If the undeformed centroidal axis is aligned along the 

/s /s 

^2 direction as shown in Figure 5.2, then £ is simply z Ji has 
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y\ 

no component in the direction. The derivative of the vector jc may 
then be written 


• ^ . 
r = {i} 


u. 


u. 


0 


(5.8) 


Remember that the vector is fixed relative to the inertial axes and 
thus the location of the undeformed centroid, z, has a zero derivative. 
The kinetic energy of the beam may now be written lin greater detail 


^ ’ k \ f ®(“? + ubdz + 4 y" Ip4»^dz 


(5.9) 


B “ B “ B 

This relationship has been simplified by eliminating terms made zero by 
the centroid definition (Equation (5.5)) and by recognizing the mass 
per unit length, m, and the mass moment of Inertia per unit length 


= /(Pl 

clA 


+ P 2 ) dm. 


(5.10) 


For the rotor, the center of mass lies on the centroid and is 
located relative to the system center of mass by the vector p . The 
derivative relative to inertial space of is The kinetic energy 

of the rotor is simplified because of terms set to zero by the defini- 
tion of the rotor center of mass. The kinetic energy of the rotor is 
written 


fp * p dm ~ 4 ML u • u -f* “ 0 )^^ • n • 

- 2 ^ -R -R 2 -R ‘-' r -R 


(5.11) 


where 
,vi 


and {~jj^ are the rotor mass and inertia dyadic respectively and 


vx 

W is the angular velocity of the rotor-^fixed axes relative to inertia: 
space. Care must be taken in the evaluation of since it is not a 


74 


first order term and since all second order terms must be retained in 
Equation (5.11). The rotor angular velocity is 


, vi vb , bi 
0) = to +0) 


(5.12) 


In this expression, the angular velocity of the rotor relative to the 
vb 

mounting frame, w , is not a first order quantity. It has a magnitude, 
Qt which is an arbitrary constant plus a first order variation and it 
has direction 
vb A 

W = fi! b ^ . (5.13) 

The angular velocity of the beam tip with respect to inertial space is 

b li. " 

jo^ . The detailed derivation of the components of is given in 

Appendix B. The second order terms must be included in the analysis. 


bi 


/ 

r * 1 

Ro 



/ 

2 



= {b}^ 1 ^ 

u' 

"l 

► + - 


V 

■ ’^R' 




(5.14) 


The components of u^^ in inertial axes are (u^ , u^ , 0) . The primes 


R, 


R. 


1 “2 

indicate spatial derivatives and the dots represent time derivatives, 


so that 


9 u. 


R. 


R, 9z 9t 


9 u 


R, 




9z 9t 


(5 .15) 


As noted above, the twist due to torsion is the angle and the sub- 
script R denotes evaluation at the rotor mounting point. The rotor 
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has a symmetric moment of inertia Ig and transverse moments of inertia 
The inertia dyadic may be written as 



(5.16) 


Keeping all second order terms, the last quantity in Equation (5.11) 
is written 


1 ,,vi 

2 - 


• Dr • + I 


2 S 


1 h{\f ^ 1 h{\f * ^ 


It is now possible to write the system kinetic energy 
T = ~ m(uj + U2)dz + ^ f dz 


1 

2 




R. 




+ i Ig 




+ 


i h{%f ^ 1 * 1 W' 


(5.18) 


The kinetic energy includes all second order terms. The subscript R 
represents evaluation of the quantity at the rotor mounting point, 
where z equals 

The strain energy of the system comes from three types of deforma- 
tion: in-plane bending, out-of-plane bending and torsion about the 

centroidal axis. The strain energy may be written 


76 


RPlPRODtrCIBILITY OF THE 

jeioo'^al page is poor 


(5.19) 



+ El 



Idz + 




where El is flexural rigidity and GJ is torsional stiffness. 

It is now possible to write the Lagrangian in terms of the scalar 
functions u^(z,t), U 2 (z,t) and \p(z,t). By taking the first variations 
of the kinetic energy in Equation (5.18) and of the strain energy in 

Equation (5.19) it is now possible to write Hamilton’s principle in the 
form. 


-jf m(u^6u^ + u^6n^)dz + ^ IpWdz 

+ M^Up(i?-p)6u^(;’^) + V2(£p)(Su2(5-p) 

- Ig h 5u^(Jl^) + Ig h Ii'(i^p)6u^(£^) 

+ Ig u|(£^)6ij;(£^) + Ig u’(S,j^)6u^(il^) 

+ nz^)6^pa^) 

~J_, V 6u" dz 


- f ^ GJ ij'Sip’ dz 


dt = 0 


(5.20) 


In order to arrive at the final formulation of Hamilton’s Principle 
it is necessary to integrate by parts. When variations of the coordi- 
nates at times t^ and t^ are set to zero, Hamilton’s Principle has the 


form 
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112 2 
£ 

2 "* 

+ [eIjU- (2j.) - M^U2a^)]6u2ai) 

+ [-Is“2<*l> - “2“2<*1> + 

+ [igRi^a^) - - Gji|i'aj)J«iK2^) 

- [eIj^u- (-«-2)]«u^(-*2)+[ei^u;i;(-^2>]s<-^2> 

-[eI 2U^' W^)] 'Su„(-Z2)+[^EIjU”(-£2)]«u^(-«2)-[GJiii' (-2^)] ailK-S^)) dt=0 

(5.21) 

Since the variations are independent and arbitrary, all quantities in 
square brackets must be set to zero. For In-plane bending, u^, the 
partial differential equation is 


+ EI^u^" = 0 

with the boundary conditions 
u”(0) = 0 
u”' (0) = 0 

EI^u”(L) = - I^li|(L) - Ig$iii,(L) 
EIiU^’ (L) = M^u^(L) 


(5.22) 


(5.23) 
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A shift of origin was made in the spatial coordinate so that -^2 became 
zero and became L, the beam length. The boundary conditions for 
in-plane bending depend on the torsion \jJ, The partial differential 
equation for torsion is 

- GJlp" = 0 (5.24) 

with the boundary conditions 

ijj* (0) = 0 

GJi|;'(L) = - I^i(L) + Igfiu^CL) (5.25) 

Again there is coupling between torsion, ip, and in-plane bending, u^, 
in the boundary conditions. This represents the gyroscopic coupling 
present in the system. Both coupling terms in the boundary conditions 
are seen to depend on the rotor spin rate, The partial differential 
equation for out-of-plane bending is 

mU2 + EI 2 U 2 "" =0 (5.26) 

with the boundary conditions 

u’^(O) = 0 

U 2 "' (0) = 0 

El2U^(L) = - IgU^(L) + IgS 

El2U^' (L) = M^U^(L) (5.27) 

The out-of-plane bending is not coupled to other system behavior, but 
the boundary conditions are time dependent. The term in Equation 
(5.27) represents a small acceleration in rotor speed which is a first 
order function of time. 
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Hamilton's Principle has been used to formulate the partial 
differential equations and the boundary conditions which govern system 
behavior. The boundary conditions represent moment and shear forces 
applied to the beam ends. The structure of the system equations demon- 
strates gyroscopic coupling between in-plane and torsional motion. 

The next section will deal with the solution of these equations by using 
separation of variables. 

5.3 Solution of the Partial Differential Equations 

In this section, the technique of separation of variables is 
employed to solve the partial differential equations which govern the 
system behavior. The first approach is to choose independent coordi- 
nates for the separation of variables. The out-of-plane motion is 
independent but the boundary conditions include time dependent terms 
and mixed spatial and time derivatives (Equation (5.27)). The separa- 
tion of variables method cannot be applied to this equation and boundary 
conditions because of their structure. Solution techniques involving 
integral transform methods such as Laplace transforms may be applicable, 
but they will be hampered by the complexity of the system and are not 
investigated here. Instead, a modified system of greater simplicity is 
introduced for out-of-plane bending. A generalized separation of vari- 
ables method is applied and the problem is transformed into one contain- 
ing a nonhomogeneous differential equation with homogeneous boundary 
conditions. This system of equations is solved using the normal mode 
method. 

Investigation of the coupled solution for in-plane bending and 
torsion shows that the choice of independent coordinates will allow 
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only the translational rigid body mode and the trivial solution to be 
admitted. The system motion is found by using a separation of variables 
with distinct spatial functions (mode shapes) but with the same time 
dependent function (modal coordinate) , The in-plane bending and torsion 
then have the same natural frequency. Solution of the partial differ- 
ential equation and boundary conditions yields "complex” mode shapes. 

This means that the modes, or eigenfunctions of the system, are func- 
tions with values which are complex numbers. The physical interpreta- 
tion of these complex mode shapes is that in-plane bending is ninety 
degrees out of phase in the time domain with torsional bending. With 
both motions harmonic, the system normal modes oscillate back and forth 
between pure in-plane bending and pure torsion. The tip of the rotor 
shaft mounted on the beam can then be seen to travel in an ellipse with 
one axis representing in-plane bending and the other axis representing 
torsion. This represents a coning of the rotor shaft about the vertical 
where the cone angle varies sinusoidally. This motion is shown in 
greater detail in the next section. 

Addressing the out-of-plane bending motion, the structure of the 
boundary conditions (Equation (5»29)) presents two degrees of complexity. 
The first property of the boundary conditions is that they involve both 
spatial and time derivatives. The presence of time derivatives normally 
causes the boundary conditions to depend on the eigenvalues of the sys- 
tem. The second property is that the boundary conditions involve the 
accelerations of the rotor and are therefore time dependent. Each of 
these complexities could be dealt with separately by known methods 
involving separation of variables, but these methods do not work for 
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the problem at hand. The methods for transforming time dependent 
boundary conditions into homogeneous boundary conditions require that 
the system does not have time derivatives in the boundary conditions. 

In order to proceed to a solution for out-of-plane bending, the physical 
system is further simplified and the time derivatives are eliminated 
from the boundary conditions. The system for out-of-plane bending is 
one where the rotor and beam, taken together, have the mass and inertia 
properties of a uniform beam. This eliminates the moment and shear 
forces applied to the beam tip by the "nonspinning" rotor. The equa- 
tion describing the system is 

mu^ + El2u'^' =0 (5.28) 

with the simplified boundary conditions 

U2"(0) =0 

U 2 "’ (0) = 0 

EIu”(L) = Igfi 

u'^’ (L) =0 (5.29) 

These boundary conditions are now free of time derivatives and may be 
transformed to homogeneous boundary conditions by 

(5.30) 


(5.31) 

The choice of S is not unique but it must satisfy the constraints 


U 2 (z,t) = w(z,t) + S(z)Igfi 

where the step function S(z) is 

S(z) = O 2 0 < z < L - E 

2 

s(z) = L - e < z < L 
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(5.32) 


S”(0) = 0 
S'” (0) = 0 
EIS''(L) = 1 
S'” (L) = 0 

The transposed partial differential equation is then 


w 


.... ^ ^ ^ 


El, 


El ®^S^ 


(5.33) 


with the homogeneous boundary conditions 


w” (0) = 0 

w'” (0) = 0 
w”(L) = 0 

w'” (L) = 0 (5.34) 

The homogeneous solution is then found using the separation of variables 
w(z,t) = g(z)y(t) (5.35) 

This yields the natural frequencies, and mode shapes, g(z) , of a 
free-free uniform beam. The nonhomogeneous solution is then found by 
modal analysis 


+ 0)^ p^(t) . N^(t) 


(5.36) 


where 


N^(t) = - H - H I_a 

n n S n S 


and 


\ X ^ 

/ 


2(EI)^*'L-e ^ 


g z dz 
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H* = r g s'"' dz 

" " El " El 

The last expression Involves the derivatives of S across the disconti- 
nuity. The resultant derivatives of Dirac "delta" functions must be 
integrated by parts, but this allows easy evaluation of the remaining 
integrals. The solution of Equation (5.36) for zero initial conditions 
is given by the convolution Integral 

n 


f 


N (t) slnw (t-T)dT 
n n 


(5.38) 


The final out-of-plane bending for our simplified model is then 
described by 


U2(z, t) 


£ g^(z)n^(t) + IgfiS(z) 


(5.39) 


The use of a step function for S raises questions concerning the 
spatial continuity of the solution and invariance of the solution to the 
location of the step. For physical reasons, the value of e would be 
kept small due to the fact that the moment is applied at the beam tip. 

The coupled in-plane bending and torsional motions are described 
by the following partial differential equations and boundary conditions 
mii^ + EI^u||^" = 0 

u5;(0) = 0 

u'l' (0) = 0 

EIu^(L) = - I^Uj^(L) - Ig^(L) 

EIu^' (L) = Mj^Uj^(L) 

(Continued) 


84 



lp\p - GJlp" =0 
ip'(O) = 0 

GJrCL) = ~ I^iCL) + Igfiuj^a) (5.40) 

If independent coordinates are chosen, then the separation of 
variables has the form 

u^(z,t) =.f(z) X(t) 

^(z,t) = h(z) V(t) (5.41) 

The first motion to be investigated is the rigid body mode correspond- 
ing to zero natural frequency. If zero initial displacements are 
assumed, then the general solutions are given by 

Ul - (Di + D2X + + D^x^)t 

^ (5.42) 


Substituting these expressions into the boundary conditions will set 
all constants to zero except D^. The rigid body motion allowed by the 
system is in-plane translation. Rotational rigid body modes are not 
permitted. 

Proceeding with the separation of variables using independent 
coordinates, the equation for the in— plane bending modal coordinate is 


X + A = 0 


(5.43) 


where 0)^ is the natural frequency. 
El, 


- — Wt f = 0 

m 1 


The equation for the mode shape is 


(5.44) 


with the boundary conditions 
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(5.45) 


f"(0) = 0 


f" (0) = 0 

El^f"(L)X = w^i^f'CDA - ioo^igfihCDv 
El^f" (L) = - 0)^ 

The boundary conditions depend on the eigenvalue, cjo^, because of 
Equation (5.43). The introduction of the Imaginary number, i, comes 
from the equation for the modal coordinate for torsion 


V + V = 0 (5.46) 

where is the natural frequency. The solution for this equation is 
io)_t 

V = ce (5.47) 

The time derivative of ^ found in the boundary conditions (Equation 
(5.40)) can then be eliminated using 

ii»(L) = ico^ h(L)v (5.48) 

By employing a similar substitution for the time derivative of the 
in-plane bending variable based on Equation (5.43), we have 


u|(L) = ico^f'CDX 

The differential equation for torsion becomes 

Ip 2 

^ GJ ^3 ^ ® 

with the boundary conditions 

h’(0) =0 

GJh' (L)v = (ja^I^h(L)v + iO)^IgP,f ’ (L)X., 


(5.49) 


(5.50) 


(5.51) 
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Let us look in detail at the boundary conditions given by Equations 
(5,45) and (5.51). It has not been possible to eliminate the time 
dependent variables from the boundary conditions. These time dependent 
variables have distinct natural frequencies. The boundary conditions 
containing the time dependent coordinates are 


[EI^f'CL) - U^Ijf'(L)]x(t) + ria)jIgSjh(L)]v(t) = 0 

jcjn'CL) - UjI^h(L)Jv(t) - CL)J\(t) - 0 (5.52) 

The only manner that these equations may be satisfied when X(t) and 
V(t) have distinct, nonzero, natural frequencies is to set the terms 
within square brackets to zero. As a consequence of this, the boundary 
conditions in Equation (5.52) become 
f '(L) =0 
f'(L) =0 
h(L) =0 

h'(L) = 0 (5.53) 

The general solutions of Equations (5.44) and (5.50) are 

f(z) = sinS^z + cos3^z + sinh3j^z + cosh3^z 

h(z) = sin32Z + Cg cos32Z (5.54) 

where 

^1 EI^ ^1 


r2 ^ _P m2 
^2 GJ “ 3 , 


(5.55) 


When these functions are substituted into the boundary conditions that 
come from Equations (5.45), (5.51) and (5.53), the solution for the con 
stants becomes 
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Cl = C2 = C3 = = Cg = = 0 (5.56) 

This represents the trivial solution corresponding to zero movement. 

In order to solve for the system behavior, a different separations 
of variables technique must be used. The modal coordinates for In-plane 
bending and torsion now are assumed to have the same natural frequency, 
0). They are separated in time by a phase delay <(>. The separation of 
variables has the new form where the coordinates are coupled 
Ui(z,t) = f(z) A(t) 

= h(z) X(t-(p) (5.57) 

The time dependent terms may now be eliminated from the boundary condi- 
tions given by Equation (5.52). The procedure is to replace V(t) by 
X(t-4>) in Equations (5.45) and (5.51). When dividing by A(t) we may 
use the relationship 

A(t) = (5.58) 

The boundairy conditions for the coupled set of equations then have the 
form 

f"(0) = 0 

f" (0) =0 
h'(0) = 0 

EIif"(L) = - a)Ignh(L)e^^‘^‘^ 

Elf'" (L) = - w^M^f(L) 

GJh’ (L) = 0)^1 h(L) + Icol fif ’ (L)e^^'^'^ (5.59) 

i. ^ 

Using the general solution forms given by Equations (5.54) and the first 
three equations given above will give the result 
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(5.60) 


This result assumes that the natural frequency is not zero. Thus the 
rigid body mode is excluded from the solutions which use Equation 
(5.60). The last three boundary conditions of Equation (5.59) deal 
with the moments and shear forces applied to the beam by the rotor. 
Using the solution forms in Equation (5.54), the boundary conditions 
are written in matrix form. 




E I j ( c o H h i * ^ L- L* o a i- ^ I.) 


EIjl-’j ( coahWjL-^cosf'j^L) 

I .j. P ^ (a 1 n IP. ;i nhi" E ) 


J<cobC«,L)« ' 


El (b inUj^L+alnhl. jL) 

+'>^(sinlijL+stnhE^r.) +u^tij,(cosl<j_L+coahEjl,) 

l/-+4 




where 


(5.61) 






(5.62) 


The characteristic equation of the system is found by setting the 
determinant of the above matrix to zero. This will yield the equation 
for the nonzero system natural frequencies. 
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sinB2L + I,j,co cos32L^ 

\/ aj ~( cosBj ^ Lcosh 3 j ^ L - l ) 

2 

cos32^L(sin3^L + sinh3j^L) 

1/2 

+Mj^(EI^m) ' C0(cos3j_Lsinh3^L-sin3^Lcosh3^L) 

0)^'^^(cos^3j^L+cos3^Lcosh32^L-sin3^Lsinh3j^L) 


- H ;3( Igfi ) COS32L 


m ( s in 3 CO s h 3 j^L+co s 3^L s inh 3 j^L ) 


+ a )^'^^( cos 3^ Lcosh 3 jL + ^^| j “ ® (5.63) 

The above equation for the natural frequency must be solved numerically. 
Once a natural frequency has been found, the corresponding mode shapes 
come from the solution of the matrix equations (5.61). The mode 
shapes are found to be 

f(z) = Cj^[sin3j^z + sinh32^z + Rj^(cos32^z+cosh3j^z) ] 


g(z) = C^R 2 e 


00832^ 


(5.64) 


where the real constants R^ and R_ are 

X Ji 


On // 

(m EIj^) (cosh32^L-cos32b) + Mj^(sin3j^L + sinh3^L) 
^ (m^EI^)^/^(sln3^L+slnh3^L) + M^(cos 33 ^L + cosh3^L) 


90 



R. = 


(cosh3^L-cos3j^)-Irj,^^^ (sin3j^L+sinh32^L)J 

(cosh3^L-cos3^L) + (sln3j^LH-sinh3 2^L) j 


|^(m^EI^)^^^(sin32L+sinh32^L) + M^\/{iij~(cos3^L+cosh32^L) (Ig^cos32L) 




^^l™(sinh3^L-sin3^L) 


/ ^ *\! 'y 

^ (cos3j^L+cosh3 


i^>] 


Ig^2 00832^ 


(5.65) 


The above mode shape for torsion is a complex mode shape. In order to 
place the proper physical interpretation on the complex mode shape, it 
is necessary to examine the displacements. 

c t \ 

^ (5.66) 

Since we know that the displacements are real, we shall take the real 
parts of the expressions in Equation (5,66) by using the relationship 

e = cosWt + i sinWt (5.67) 

The displasements then become 

u^ = C^[sin3^z+sinh3j^z+R^(cos3j^z+cosh3^z) ]cosO)t 

ip - C^[R2Cos322]sinO)t (5.68) 

The arbitrary time lag, (}>, drops out of the equations. The equations 
for the displacement show that the in-plane bending is always ninety 
degrees out of phase with the torsional bending. This follows directly 
from the phase relationship of the sine and cosine functions. It is 
this phase relationship between in-plane bending and torsion which 
shows the system to have gyroscopic coupling. The system behavior may 
now be calculated by solving the frequency Equation (5.63) and then 
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calculating the mode shapes by (Equation (5.64)). This work must be 
accomplished numerically due to the complicated structure of the equa- 
tions. 


5.4 Numerical Solutions 

In this section, we shall formulate an example problem and solve 
for the coupled in-plane and torsional motion. In this manner, we shall 
be able to examine the gyroscopic coupling exhibited by this system for 
a specific case. This should Increase our understanding of the dynamic 
interactions involved for this class of systems. 

We start by choosing physical parameters for the system. The 
beam will be 10 meters long with an in-plane bending fundamental fre- 
quency equal to 1 Hertz and a torsional fundamental frequency equal to 
5 Hertz. The mass of the beam will be 150 kg. and it has a cross 
section which is one meter square. The relationships for the funda- 
mental frequencies of a free-free uniforms beam are 



These expressions are inverted in order to solve for the values of 
flexural rigidity and torsional stiffness 

El = 11820. 

sec 
, 3 

GJ = 25000. (5.70) 

sec 

The value for mass moment of inertia per unit length, 1^, is 2.5 kg-m. 
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The rotor has a mass of 5.236 kg, a symmetric moment of inertia of 

2 2 
.2909 kg-m and a transverse moment of inertia of .5818 kg-m . The 

rotor spin speed is fixed at values between 0 and 1000 radians per 

second. This value is varied parametrically. 

The first result to be shown is the dependence of the system 
natural frequencies on the rotor spin speed. Solutions of the frequency 
equation (Equation (5.63)) were obtained numerically for rotor speeds 
between 0 and 1000 radians per second. The first two natural frequen- 
cies are plotted in Figures 5.3 and 5.4. In order to examine the mode 
shapes themsfilves, the rotor spin speed was fixed at 1000 radians per 
second and the mode shapes were calculated using Equation (5.64), 

These mode shapes are shown for the first five natural frequencies in 
Figures 5.5 to 5.9. The first mode (Figure 5.5) involves much more 
rotation than bending. It is related to the rotational rigid body mode 
that would result if the beam were rigid. This is confirmed by noting 
in Figure 5.3 that for zero spin speed this mode has zero natural fre- 
quency. The last aspect of the system behavior to be examined involves 
the rotor motion. For small in-plane displacements, the rotor shaft 
will rotate an angle equal to the slope of the in-plane mode shape 
evaluated at the beam tip. The tip of the rotor shaft will then rotate 
in an elipse as shown in Figure 5.10. The rotor is spinning counter 
clockwise, but the rotor tip can cone in either direction. The direc- 
tion of the coning motion is determined by the moment applied to the 
rotor (negative of the moment applied to the beam) and is thus a function 
of the curvature (second derivative) of the mode shape evaluated at the 
beam tip . 
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Figure 5.7. Third Mode Shapes. 
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Figure 5.9. Fifth Mode Shapes. 
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CHAPTER 6 


NUMERICAL RESULTS FOR THE FINITE ELEMENT MODEL 
6.1 Chapter Summary 

This chapter presents a solution using numerical integration for 
the finite element model equations of motion. The physical system 
chosen for the example problem consists of a uniform beam with a con- 
stant speed, axisymmetric rotor mounted at the tip. This is the same 
system analyzed in the previous chapter which used a continuum approach 
that required the solution of the governing partial differential equa- 
tions. The verification strategy consists of a comparison of the 
results of the numerical integration of the finite element equations of 
motion with the closed form solution attained by the continuum analysis. 
The Initial conditions are calculated using the closed form solution 
and correspond to a system normal mode for the coupled in-plane bending 
and torsional motion. 

The modal analysis is performed on the system with the rotor 
"frozen" (zero spin rate). In the construction of the mass and stiffness 
matrices, constraint relationships arise since each mass element of the 
beam has only three degrees of freedom. The constraint relationships 
are accommodated through the Introduction of displacement compatibility 
matrices when the final matrices are assembled from their component 
parts. The stiffness matrix may also be calculated from the continuum 
strain energy expression by the adoption of finite differences. An 
eigenvalue analysis of the final mass and stiffness matrices produces 
a set of natural frequencies and real modes (eigenvectors). 
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The numerical Integration is accomplished by a fourth-order Runge- 
Kutta algorithm. The only variables Integrated are the Tisserand frame 
variables and the modal coordinates corresponding to the lowest nonzero 
frequencies for both in-plane bending and torsion. This represents a 
truncation of the number of coordinates from thirty to four (no out-of- 
plane motion is initiated) . The results of the integration show good 
agreement with the closed form solutions. The chapter closes with a 
discussion of the truncation procedure. Emphasis is placed on the inter- 
relation of the coordinates that results from gyroscopic coupling. 

Modes which are strongly coupled should be treated as a unit in the 
truncation procedure; that is, they should be kept as a group or 
truncated as a group. 

6.2 Example Problem and Modal Analysis 

This section will describe the physical system to be modeled and 
the constraint relationships involved. It will then construct the mass 
and stiffness matrices to be used in the modal analysis. The physical 
system is the same as the one described in Chapter 5. The difference 
here will be the method by which a mathematical model of the system is 
created. For the previous analysis, the system was modeled as a uniform 
continuum and partial differential equations with boundary conditions 
resulted. Ordinary differential equations were obtained by introducing 
a separation of variables. For the finite element model, ordinary 
differential equations are obtained directly by breaking the system into 
a number of discrete, finite elements. 

To review the selection of parameters, the beam is 10 meters in 
length, 150 kg in mass and has a 1 meter square cross section. The 


104 



3 2 

flexural rigidity, El, is equal to 11820 kg-m /sec and the torsional 

3 2 

stiffness, GJ, is 25000 kg-m /sec . The rotor has a mass of 5.236 kg, 

2 

a symmetric moment of inertia of .2909 kg-ia and a transverse moment 

2 

of inertia of .5818 kg-m . The rotor spin speed is fixed at 10 radians 
per second. 

The beam is divided into fifteen equal rigid bodies that are 

connected by massless beams as shown in Figure 6.1. Each rigid body 

has a mass of 10 kg and a moment of inertia about the third axis of 
2 

1.6667 kg-m . The other moments of inertia are set to zero in order to 

eliminate rotary inertia effects. The system moments of inertia are 

2 

1262.8, 1263.1 and 25.6 kg-m about axes one, two and three respectively. 
Each massless beam has a length of two thirds of a meter. 

We shall be concerned with the calculation of natural frequencies 
and mode shapes for in-plane bending and torsion. The plane mentioned 
here contains the rotor symmetric axis and the undeformed beam centroid. 
Torsion is about the beam centroid. If we consider in-plane motion 
first, we see that each rigid body will have three degrees of freedom; 
transverse displacement and longitudinal displacement of the center of 
mass and rotation about axis number two. The first constraint simply 
sets the longitudinal displacement to zero by assuming that the beam is 
liiextensible. The second constraint relates the rotation to the trans- 
verse displacement. This is the same curvature constraint which was 
used in the continuum analysis (Equation (5.14)) and is derived in 
Appendix B. 



( 6 . 1 ) 
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If this constraint is applied at the midpoint between rigid bodies and 
if the angle is taken to be the average value for the two neighboring 
bodies, the constraint relationships take the form 




,i+l 


H-1 


“ u. 


Az 


i = 1 


,n-l 


( 6 . 2 ) 


where Az is the distance between body centers of mass and is here equal 
to two thirds of a meter. This relationship can be solved for and 

yields 

' “l + r ef’) n-l (6.3) 

This same result may be attained using a different perspective. Let 
the bodies be connected at the midpoints between them by a line hinge, 
then the requirement for compatible displacements at the hinge point 
will yield 

“i 2 ^1 ~/^i+l ” 2 ^1+1 ^ “ l,...,n-l (6.4) 

where the rotations are assumed to be small. This will produce the same 

result (Equation (6.3)) as the curvature constraint. There are now 

(n-1) constraint relationships and they may be used to eliminate the 
2 ' 

variables u^ through u^. In matrix form, the constraint relationships 
may be expressed as 
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In compact form, this may be rewritten I 

I 

{q> = [C]{q'} (6.6) I 

where [C] is referred to as a compatibility matrix. The mass matrix d | 

£ 

will be found by applying a transformation involving the compatibility I 

matrix to the general mass matrix involving all displacements and rota— | 

tions. The general mass matrix may be written in the form ! 

■ ' ' . t. 
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where guidance has been taken from Equation (4.38) and m^ is the mass 
of body i and is the rotor transverse moment of inertia. The 
inertias of each body about axis number two have been set to zero in 
order to eliminate the rotary inertia effect. The final mass matrix 
has the form 

[MV] = [C]^[M][C] (6.8) 

where [M* ] is symmetric but is no longer diagonal. The stiffness 
matrix for in-plane bending is constructed from an assemblage of simple 
component matrices by using compatibility matrices. The simple compo- 
nent matrix is constructed for each massless beam by using assumed mode 
shapes. Each member of this matrix is found by evaluating the expres- 
sion 

ri ■ 

EI(})^4>^' dx (6.9) 

The assumed mode shapes, (})^, chosen here arise from a general third 
order- polynomial that is forced in turn to represent unit displacements 
and rotations about each end of the beam. For unit displacement about 
the end at x = 0, the mode shape is 

- 2(^)' (6.X0) 

For unit rotation about this end, the mode shape is 

<|>2(x) = (6.11) 

For unit displacement about the end at x = £, the mode shape is 

*3 ’ ( 6 . 12 ) 
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For unit rotations about this end, the mode shape is 

.2 




(6.13) 

Using Equation (6.9) will yield the component stiffness matrix for each 
massless beam element 




12 

6Z 

-12 

6Z 

61 


-6£ 

2Z 

-12 

-6Z 

12 

-6Z 


2 

-6Z 

AZ 


(6.14) 


The general stiffness matrix is then written in the form 


[K] - 


[kiJ 


[k2l 




[k ,] 

n-1 


(6.15) 


An additional compatibility matrix is also Introduced to set equal the 
displacement and rotation at x = il of beam i to the displacement and 
rotation at x = 0 of beam i +1, respectively. 


[C’] = 
(4(n^l)x2n) 


1 

0 

0 

0 

0 

0 

0 

0 


0 

1 

0 

0 

0 

0 

0 

0 


0 

0 

1 

0 

1 

0 

0 

0 


0 - - - 

0 

0 - - - 
1 0 - 


0 

1 

0 

0 


0 

0 

1 

0 


0 

0 

0 

1 


0 - - - ~ 


0 0 0 - 


(6.16) 



Note that the structure of [C] consists of (n-1) four by four identity 
matrices with each one overlapping the last two columns of its predeces- 
sor, The final stiffness matrix may then be written in the form 

[K'] = ([C’][C])'^[K]tC’][C] (6.17) 

A different perspective may be taken in order to derive the strain 
matrix. To do this, we start with the continuum expression for the 
strain energy 

•L /a^uA 

(6.18) 

where L is the length of the entire beam. Introducing the curvature 
constraint (Equation (6.1)) yields 

2 



« = 2 i 


(6.19) 


This expression may now be evaluated using a finite difference approach, 
First the substitution is made 


3^2 ^ 


( 6 . 20 ) 


8 z Az 

where Az is the distance between mass centers. The strain energy may 
then be written 


u a ^ 


2 2 
i=l (Az) 


(p?" - 32 J 


(6 . 21) 


Expressed in matrix form, the stiffness matrix, [K^], may be introduced 

V = (6V22) 

where { 62 ) is the column matrix of rotations. In detail, the stiffness 


matrix has the structure 
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1 -1 0 0 
-1 2-10 


= El 
(Az)‘ 


0 -1 2 -1 

0 0 -1 1 


(6.23) 


An interesting result is that the strain energy does not involve 
the displacement This means that the matrix [K’ ] found in Equa- 

tion (6.17) must be singular; that is, it must possess a zero first 
row and column. The matrix multiplications indicated in Equation 
(6.17) were performed numerically and the stiffness matrix was found 
to be singular. The first row and column were found to be zero and 
the populated submatrix was equal to [K^] found in Equation (6.23). 

This zero row singularity will yield a zero eigenvalue which will cor- 
respond to the translational rigid body mode. Since the rank of [K'j 
is degenerate by 2 (two less than order of the matrix) another zero 
eigenvalue will be present and corresponds to the rotational rigid body 
mode . 


The free-free modal analysis of the in-plane bending is performed 
by solving the eigenvalue problem 


([K'] “ W^[M']){q’} = {0} (6.24) 

The resulting eigenvectors {q|} for distinct natural frequencies 
are orthogonal since the mass and stiffness matrices are symmetric. 

{qi}'^[M']{qp = 0 ± ^ j 

0)^ oJj (6.25) 

If the eigenvectors are arranged by columns in a modal matrix [<j) ’ ] , and 
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if the eigenvectors corresponding to the same natural frequency are 
made orthogonal by choosing their proper linear combinations, then the 
system may be normalized to obey the relationship 

= [E] (6.26) 

where [E] is the identity matrix. The constraint relationships may be 
reapplied to the resultant mode shapes using Equation (6.6) and the 
resulting mode shapes will represent compatible displacement and rota- 
tion. In terms of the modal matrix, we express this operation as 

m = [C][(|)*l (6.27) 

These resulting modes may now be input into the finite element model 
equations of motion. Plots of the first three nonzero frequency modes 
are shovm in Figure 6.2. 

The modal analysis of the torsional motion is much simpler due to 
the lack of constraints. The mass matrix is diagonal and has the form 


(6.28) 

where is the moment of Inertia of body 1 abo-.t axis number three. 

The stiffness matrix is calculated from the continuum strain energy. 


[M'3 = 


■3 s 


3 



(6.29) 


by the introduction of the finite difference 


3fi, 


oz 


. 1+1 .1 

^3 ~ ^3 

Az 


X 1 , . . . , n— 1 I 


(6.30) 
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The strain energy may then be written 

u = i E (ef ' - iif 

1=1 (hz) \ V 

If we write the strain energy in matrix form, we introduce the 
stiffness matrix, [K'J 


(6.31) 


U = \ {B3}’^[k']{33} 


(6.32) 


where { 33 } is the column matrix of rotations about axis three for all 
bodies. The structure of [K*] is 


[KM .= 


GJ 


(Az) 


2 .' 


1 

-1 

0 

0 


-1 

2 


0 0 

-1 0 

» 

= - 0 

~ - 0 


-1 

0 


0 

- - 0 

2 -1 

-1 1 


(6.33) 


It is now possible to proceed with the free-free modal analysis of the 
torsional motion by solving the eigenvalue problem 


([KM - W^[MM){33} = {0} 


(6.34) 


[4>]'^[m’][4>] = [E] 


The modal matrix, [(})], is formed by the eigenvectors arranged in columns 
and is normalized to obey the relation 

(6.35) 

The first three nonzero modes are shown in Figure 6.3. One rigid body 
mode corresponding to uniform rotation about axis three is present. 

The modal analyses is now complete. The real mode shapes and 
natural frequencies for in-plane bending and torsion found here may be 


115 



f = 4.880 Hz 





< f = 14.471 Hz 


Figure 6.3. Torsion Mode Shapes. 
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used In the finite element model equations of motion. The next step is 
the actual integration of the equations of motion. 

6.3 Numerical Integration and Truncation Procedure 

This section presents the results of a numerical integration of 
the finite element model equations of motion for a system with a con- 
stant spin rate. The initial conditions chosen for the integration 
correspond to a system normal mode and should therefore produce periodic 
motion. The initial conditions (Figure 6.4) were calculated using the 
closed form solutions resulting from the continuum analysis of the 
previous chapter. 

The initial modal coordinates and velocities are calculated using 
the inverse modal matrix. This may be represented in a very simple 
form due to the normalized orthogonality condition (Equation (4.47)) 

(W]^[M]) [<})] - [E] (6.36) 

This then requires that the first two matrices in parentheses be equal 
to the inverse modal matrix 

[M] (6.37) 

The initial modal coordinates are then 

The resulting initial modal coordinates and velocities corresponding to 
nonzero frequency modes are shown in Table 6.1. The Initial Tisserand 
frame variables are simply the weighted mean values of the initial 
angular variables. The Initial value is -.17515 E-4, the initial §2 
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TABLE 6.1 


INITIAL MODAL COORDINATES AND VELOCITIES 
Torsion In Plane Bending 


Natural 
Erequency (Hz) 

Modal 

Coordinate 

Natural 
Frequency (Hz) 

Modal 

Velocity 

4.880 

.47711837 

E-5 

0.964 

-.72829664 E-1 

9.718 

-.11208203 

E-5 

2.713 

-.36030891 E-3 

14.471 

-.47018966 

E-6 

5.427 

-.15179138 E-3 

19.091 

-.24678593 

E-6 

9.165 

.39719278 E-4 

23.527 

-.14601846 

E-6 

14,023 

.32309908 E-4 

27.727 

-.92882146 

E-7 

20.138 

.78941230 E-5 

31.641 

-.62197614 

E-7 

27.741 

.14494639 E-4 

35.218 

.43088221 

E-7 

37.273 

-.60349703 E-5 

38.415 

.30428964 E-7 

49.649 

.33760443 E-5 

41.190 

-.21900632 

E-7 

66.720 

-.56705903 E-5 

43.511 

-.15267688 

E-7 

92.398 

.44775661 E-5 

45.346 

.10585609 

E-7 

136.712 

.23686443 E-5 

46.673 

-.64728738 

E-8 

236.365 

-.69689122 E-7 

47.477 

-.31017748 E-8 

723.839 

-.13253566 E-5 
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value is -.17337 E-3, the remaining Tisserand frame values are initially 
zero . 

An examination of the initial modal coordinates and velocities 
will show that the mode shape represented in Figure 6.4 is dominated 
by the first In-plane bending mode (.964 Hz) and by the Tisserand frame 
variable, 9^, representing uniform rotation about axis number three 
(corresponding to the torsional rigid body modes). The first torsional 
mode (4.880 Hz) is also used in order to provide a variation of torsional 
rotation along the length of the beam. The numerical simulation will be 
performed for the Tisserand frame variables and the lowest frequency 
bending and torsion modes. This represents a truncation of the remain- 
ing twenty six coordinates. 

The integration was performed by a fourth-order Runge-Kutta 
algorithm with a step size of .01 seconds. The integration time span 
was set equal to 1.06 seconds, the period of the normal mode used for 
the initial conditions. 

The first aspect of the results of the numerical integration to be 
examined is the simulation of the natural frequency of the normal mode. 
The times corresponding to the quarter cycle points are shown in Table 
6.2. For the closed form solution, the bending and torsion have equal 
natural frequencies, but for the numerical integration, u^ and 3^ have 
unequal natural frequencies. The in-plane bending is well behaved and 
has a period of 1.04 seconds. The torsional motion, follows the 
periodic motion less precisely. A comparison of the displacements is 
shown for the quarter cycle points in Figures 6.5 to 6.8. Again, the 
in-plane bending shows very excellent agreement, while the agreement for 
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TABLE 6.2 


COMPARISON OF NATURAL FREQUENCIES 


t (degrees) 

Times of 

Occurrence 

(seconds) 


Closed Foinn 

Numerical Integration 


Solution 


^3 

0 

0.0 

0.0 

0.0 

90 

0.27 

0.26 

0.27 

180 

0.53 

0.52 

0.50 

270 

0.80 

0.78 

0.72 

360 

1.06 

1.04 

> 1.06 


torsion is not as precise. It is important to note, however, that the 
magnitude of the torsional motion is very small, and the imprecision 
does not result in a great error in the position of each finite element. 
The last aspect of the numerical integration to be examined is the 
motion of the tip of the rotor shaft. As shown in Figure 6.9, the rotor 
shaft follows the correct path as determined from the closed form solu- 
tion with only a small imprecision. 

The agreement of the numerical integration results with the closed 
form solution is good. The deviations which do exist are small and 
may be attributed to three sources. First there is the difference in 
the mathematical models used; a continuum analysis should agree with a 
finite element analysis only in the limit. For a small number of bodies 
as chosen here, a difference in results will occur. The second source 
is the truncation procedure. The number of coordinates was reduced from 
thirty to four. The truncation alters the angular momentum of the 
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“• CLOSED FORM SOLUTION AT .27 SEC 
O NUMERICAL SOLUTION AT .26 SEC 


Figure 6.5, Numerical Results at Quarter Period Point 
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Figure 6.8. Numerical Results at One Period Point. ' | 
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system, which probably causes the shift of the center of the locus of 
the tip of the rotor shaft seen in Figure 6.9. The computational ease 
of Integrating the truncated set of finite element equations of motion 
is one of the main advantages of the formulation. In fact, if the 
number of modes is increased, not only do the computations become more 
numerous, but they also experience numerical difficulties arising from 
small differences of large numbers. The third source is a small dif- 
ference in the mounting location of the rotor shaft. For the continuum 
analysis, the rotor is mounted at the tip, at L = 10 m. In order to 
collocate the centers of mass of the rotor and the last rigid body, the 
rotor shaft is mounted at L = 9.333 m, in the finite element model. 

For this numerical example, the truncation was performed in an 

ad hoc manner by examining the initial modal coordinates. A more formal 

procedure would necessitate the examination of the modal gyroscopic 
T 

coupling matrix, ^ G(p\ The. structure of the coupling matrix would 
indicate which modes were strongly coupled for the case at hand. 

Strongly coupled modes correspond to components of a complex mode shape. 
They must be truncated or incorporated in the system of equations as a 
unit. Partial truncation of strongly coupled modes would introduce 
unnecessary truncation errors. A further point concerns the analysis 
of the sensitivity of vibrations to excitation sources. True natural 
frequencies of the system are to be avoided or somehow compensated for 
by the control system. It is Important to realize that the natural 
frequencies present in the equations of motion are not the true system 
natural frequencies when the rotor spin rates are nonzero. For nearly 
constant spin rates, the true natural frequencies correspond to the 
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complex mode shapes. If the spin rates vary greatly, a range of natural 
frequencies would have to be included in the analysis of the effects of 
excitation sources. An example of this behavior was given in the 
previous chapter in Figures 5.3 and 5.4. The relationship between 
complex modes and real modes must then be made by identifying the domi- 
nant terms when the components of a complex mode are expressed as linear 
combinations of real modes. It is obvious that the truncation procedure 
for the finite element model equations of motion will be complicated for 
true spacecraft systems. The above must certainly be expanded upon by 
further investigation. 

This completes the numerical investigation of the finite element 
equations of motion. The principal conclusion is that it is possible to 
represent the behavior of a complex normal mode by a set of real modes 
in a coupled set of equations of motion. 
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CHAPTER 7 

SUMMARY, CONTRIBUTIONS AND FUTURE WORK 
This chapter will tie together the results of the previous work 
and present conclusions* This dissertation attempted to deal with the 
problems presented by the next generation of spacecraft. We are enter- 
ing a new period in the analysis of the effects of spacecraft flexibil- 
ity. Until now, the effects of flexibility for real spacecraft systems 
have consisted mainly of those produced by flexible appendages attached 
to central rigid bodies. This research addresses the more general case 
where the structural flexibility is distributed throughout the body. 

The problem is significantly complicated by the presence in the system 
of momentum exchange controllers. 

The first part of the dissertation consists of a detailed 
investigation of floating reference frames. These frames are defined 
in such a way that they somehow follow the overall motion of the body. 
When deformations are ejcpanded relative to these frames, they are small 
This makes a first order vibration analysis possible and thereby yields 
a significant advantage. One possible drawback to their use is the 
introduction of the constraint relationships which define their motion. 
An investigation is made of the constraint relationship which defines a 
floating frame by setting to zero the linear and angular momentum 
relative to the frame. This is the Tisserand frame, named after a 
nineteenth century French astronomer. A special case of this frame is 
found which uses the free-free modes of the unconstrained system, and 
which satisfies the constraint relationship if the rigid body modal 
coordinates are simply set to zero. Since the remaining modal 
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coordinates are independent due to properties of their mode shapes, it 
is called the Mode Shape Constraint. It is possible to expand these 
ideas to cover a deformable body containing spinning rigid rotors by 
defining a Modified Tisserand Constraint. 

The main body of the dissertation deals with the formulation of 
equations of motion for spacecraft with distributed flexibility and 
momentum exchange controllers. We are interested in addressing the 
coupled structure-actuator dynamics. The control variables for these 
systems would be the rotor spin speeds. For the specific application 
studied, the rotor torques must formally be treated as small, first 
order quantities. For engineering systems. It may be possible to relax 
this requirement, and the rotor spin speeds will then vary greatly. If 
complex mode shapes, for which the eigenvectors are complex numbers, 
were used in the formal analysis, the resulting equations of motion 
would not include the spin speeds of the rotor as explicit functions. 

Any relaxation of the first-order torque requirement would probably 
require the staging of complex modes for different ranges of spin speeds 

The technique developed here uses the real modes of the uncon- 
strained system where the rotors have been "frozen", that is, they have 
zero spin rates . The resulting derivation requires that the modal 
coordinates are coupled, but it is possible to truncate the system to a 
smaller number of coordinates and still achieve accurate results. 

A main advantage of the final set of equations of motion is that 
the control variables, the rotor spin speeds, appear as explicit terms 
in the equations. The Modified Tisserand Constraint is used in the 
formulation in order to provide an example of its properties. The 
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Modified Tisserand frame may also make it easier to relax the constraint 
that rotor torques are small, first order quantities. 

Considerable attention is paid to developing physical insight into 
system behavior and the numerical verification of results. A continuiam 
analysis is used to formulate closed form solutions for a simple example 
system consisting of a uniform beam with a nearly constant speed, 
axysymmetric rotor mounted at the tip. The governing partial differen- 
tial equations and boundary conditions are obtained in a unified, sys- 
tematic fashion by the application of Hamilton's Principle. Closed 
form solutions are obtained through the use of generalized techniques 
for the separation of variables. For a specific choice of parameters 
and a constant spin speed, system normal modes are c^slculated. The 
rotor shaft is found to cone with a sinusoidally varying cone angle. 

The numerical verification technique of the finite element equa- 
tions of motion also uses this simple physical system. After the real 
mode shapes have been calculated for the system with a "frozen” rotor, 
the initial conditions are chosen with the aid of the closed form solu- 
tions so that they produce a periodic system normal mode. The truncated 
equations of motion are integrated numerically for one period and the 
results compare well with those predicted by the continuum analysis. 

The research undertaken in the dissertation has contributed to the 
understanding and solution of an engineering problem of current interest 
Specifically, the contributions of the dissertation are: 

1. Formulated a set of finite element equations of motion suitable 
for the analysis of spacecraft with distributed flexibility and 
momentum exchange controllers . A significant advantage of this 
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formulation is that the control parameters (rotor spin speeds) 
appear explicitly in the equations, and that the coordinates 
can be truncated. 

2. Gained insight into the physical behavior of this class of 
dynamical system by obtaining closed form solutions using a 
continuum analysis. 

3. Expanded the understanding and usage of the floating Tlsserand 
reference frame. Made it possible to bring the advantages of 

a rigorous first order vibration analysis to a system experienc- 
ing large rotations. 

It seems appropriate that this dissertation should close with an 
enumeration of work still left to be done. It seems to be in the nature 
of things that each answer only necessitates further questions. The 
areas of future research connected with this dissertation are: 

1. Control system analysis of the finite element equations of 
motion. This is a rich field for future work and involves 
the assessment of questions of stability, observability and 
controllability. Specifically, questions of where to put 
sensors and actuators , and how many are needed, will have to 
be answered. Transfer functions should also be developed. 

2. A method of synthetic modes should be developed in order to 
compensate for the error in system angular momentum caused 
by truncation. 

3. The finite element equations of motion could be formally 
extended to cover large rotor torques and large rotations 
of the Tlsserand frame. This could aid in assessing the 
possibility of relaxing the first order torque requirement 


132 



lor the simpler system of equations presented in Chapter 
four. 

Further analysis of the truncation procedure. Formal 
procedures need to be developed by expanding on the ideas 
presented here. Sensitivity of the model accuracy to trun- 
cation needs to be assessed. 

The class of dynamical systems covered could be expanded. 
One extension would involve the addition of an articulated 
member to the system. 
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APPENDIX A 


NOMENCLATURE 


Ref « Eg , 


Constraint relationship coefficients for 
a Pfaffian constraint 

Body i in the finite element model 

Vector basis fixed in the mass element locally 
attached frame. Commonly at the rotor mounting 
point 

Compatibility matrix 

Constant 

Constant 

Beam cross section normal to centroid 
Beam flexural rigidity 
Forcing matrix 

Spatial variable for in-plane bending 

Vector basis fixed in the system floating frame 
commonly the Tisserand frame 

Bending fundamental frequency 

Torsion fundamental frequency 

Gyroscopic matrix, skew symmetric 

Beam torsional stiffness 

Spatial variable for out-of-plane bending 
System angular momentum 
Rotor angular momentum 
Spatial variable for torsion 

System inertia dyadic about system center of mass 


( 2 . 18 ) 

( 4 . 7 ) 


( 6 . 6 ) 

( 5 . 54 ) 

( 5 . 42 ) 

( 5 . 5 ) 

( 5 . 19 ) 

( 4 . 37 ) 

( 5 . 41 ) 

( 5 . 69 ) 

( 5 . 69 ) 

( 4 . 37 ) 

( 5 . 19 ) 

< 5 . 35 ) 

( 2 . 1 ) 

( 4 . 17 ) 

( 5 . 41 ) 

( 2 . 3 ) 




[I] 


[I^] 

[iR] 

[K] 

L 

{L} 


jU 

[M] 



Ref . Eq . 


Deformable body inertia dyadic about the system 
center of mass 

(3.32) 


Rotor inertia dyadic relative to rotor center 
of mass 

(3.33) 


Rotor inertia dyadic about the system center 
of mass 

(3.33) 


Polar mass moment of inertia per unit length length 

(5.9) 


Rotor symmetric axis principal moment of inertia 

(5.16) 


Rotor transverse axis principal moment of inertia 

(5.16) 

22*^33 

Moments of Inertia of the undeformed system 

(3.22) 

13’^23| 

T 1 

sProducts of inertia of the undeformed system 

(3.22) 

31* 32> 

1 

System inertia matrix 

(4.34) 


Inertia matrix of body i 

(4.14) 


Inertia matrix of body to which the rotor is 
attached 

(4.34) 


Inertia matrix of the rotor 

(4.14) 


Vector basis fixed in inertial reference system 



Stiffness matrix-symmetric 

(4.36) 


Lagrangian of the system 

(5.1) 


Beam length 

(5.23) 


Lagrange multiplier coefficient matrix 

(4.37) 


Distance between beam free end and system 
center of mass 

(5.20) 


Distance between rotor mounting point 
and system center of mass 

(5.20) 


System mass 

(3.14) 


Mass matrix-symmetric 

(4.37) 
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Rotor mass (3.33) 

Mass per unit length for a beam (5.9) 

Location of a beam mass element in the plane 

normal to the centroid (5.4) 

Scalar components of £ in the beams locally 

attached coordinate system with the origin 

at the centroid (5.6) 

Matrix of generalized forces (4.36) 

Generalized force corresponding to generalized 
coordinate k. (4.26) 

Matrix of displacement and rotation coordinates (4.35) 

Discrete coordinate coupling matrix (4.37) 

Constants (5.68) 

Location of the beam centroid relative to 

the system center of mass (5.4) 

Location of the deformable body center of mass 

relative to the system center of mass (3.44) 

Location of generic mass element relative to 

body i center of mass (4.5) 

Spatial step function (5.30) 

Components of w the coordinate 

system of body R (4.14) 

Kinetic Energy (2.5) 

System external torques (4.30) 

System potential energy (5.2) 

Deformation of a generic mass element. 

For a beam, the deformation of the centroid. (2.10) 

Displacement of body 1 center of mass 

due to deformations (4.5) 


, {u^} 


Ref . Eq . 

Components of u in the floating frame coordinate 
system 

(4.14) 


Deformation of the rotor mounting point 

(5.11) 


0) Components of in inertial coordinate system 

(5.14) 

(ui.U2>“3 

) Components of u in the floating frame coordinate 
system 

(2.13) 


{v} 

Vector bases fixed in the rotor frame 


w 

Component of out of plane bending 

(5.30) 

(x,y,z) 

Components of £ in the floating frame 



coordinate system 

(2.13) 

{x^} 

Components of £^ in floating frame coordinate 
system 

(4.23) 

{3^} 

i 

Components of ^ ^ in body i coordinate system 

(4.14) 

%’^2 

Constants related to angular velocity 

(5.55) 

[A] 

Modal coordinate coupling matrix 

(4.55) 

e 

Small number used in step function 

(5.31) 

[?] 

Diagonal matrix of modal percent of critical 



damping terms 

(4.54) 

{n} 

Matrix of modal coordinates 

(4.48) 


Modal coordinate 

(2.14) 


Rigid body modal coordinate 

(3.13) 

{6} 

Tisserand frame variables 

(4.13) 

A 

Time variable for in-plane bending 



(later for torsion also) 

(5.41) 


Lagrange multiplier term corresponding to 


{A} 

constraint relationship s. 

(4.26) 

Lagrange multiplier matrix 

(4.37) 

y 

Time variable for out of plane bending 

(5.35) 
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V 


Time variable for torsion 


Ref. Eg . 
(5.41) 



£ 

Location of the generic mass element relative 
to the center of mass 

(2.1) 

£ 

Time derivative of £ relative to an inertial 
reference frame 

(2.1) 

O 

£ 

Time derivative of £ relative to floating 
reference frame f. 

(2.2) 

£ 

Location of the generic mass element in the 
undeformed state. Fixed in the floating 
reference frame 

(2.10) 


Location for the undeformed state of the generic 
mass element of the deformable body relative to 
the deformable body center of mass . Fixed in the 
floating reference frame 

(3.30) 


Undeformed location of body i center of mass 
relative to the system center of mass. Fixed 
in the floating reference frame 

(4.5) 

■^R 

Location of the rotor center of mass relative 
to the system center of mass 

(3.33) 

(P2>P2*P2^ 

Components of £ in the floating frame coordinate 
system 

(2.7) 

m 

Rectangular modal matrix excluding rigid body modes 

(4.47) 


Mode shape 

(2*15) 

if 

Rigid body mode shape 

(3.7) 


Components of In the floating frame coordinate 

system 

(2.14) 


Deformational mode shape corresponding to a 
nonzero natural frequency 

(3.10) 

'I) 

Rotation about the beam centroid due to torsion 

(5.7) 


Rotation of the rotor about the beam centroid 
due to torsion 

(5.14) 

a 

Rotor spin speed 

(5.13) 



r 
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Ref . Eg . 



System natural frequency 

(5.36) 

M 

bl 

0) 

Diagonal angular velocity squared matrix* 
Nonzero frequencies only 

(4.47) 

Angular velocity of beam locally attached frame 
relative to inertial frame 

(5.12) 

bf 

0) 

Angular velocity of the beam locally attached 



frame (commonly at the rotor mounting point) 
relative to the system floating frame 

(3.36) 

b^f 

Angular velocity of body i frame relative to 


w 

(4.9) 

Tisserand frame 

fi 

Angular velocity of the floating frame 



(2.2) 

relative to an inertial frame 

vb 

03 

vbR 

i£ 

Angular velocity of rotor fixed frame 
relative to the beam locally attached 
frame at the mounting point 

(3.36) 

Angular velocity of the rotor fixed frame 
relative to the body R frame 

(4.10) 


Angular velocity of the rotor frame relative 
to the floating frame 

(3.33) 




APPENDIX B 


COORDINATE TRANSFORMATIONS 
1-2-3 COORDINATE TRANSFORMATION 


{b} = 


COS0^ 

sin0 

-sinOg 

COS0 

0 

0 


Of 


COS02 0 -sin02 


_sln02 0 COS02J 


1 0 
0 COS0 
0 -sin0 


1 

^ COS0^J 


{ 1 } 


{b} 


COS02COS02 


cos0j^sin0 2 


sin0^sin02 


+sin0 j^sin0 2 COS 0 ^ -cos0 ^sin0 2 cos 0 ^ 


■cos02sln02 cos0^cos02 


sin0^cos02 


{ 1 } 


-sin0^sin02Sin02 +cos0^sin02sln02 


s±n0, 


-sin0^cos02 


bi • ^ 

w = (0^cos02cos0^ + 02Sin02) 

+ (-5^cos02sln02 + 02COS02) ^2 

+ ( 0 ^sin 02 + 63 ) b^ 


4- (02cos0j^ - 02Sin0^cos02) ±2 
+ (02sin0^ + 02cos0^cos02) I.3 


First Order Approximation 


cos0^cos02 


{b} = 


1 

-0, 

0 . 


0, 

1 

-0. 


{1} = ([E]-[e]){i} 
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{ e, \ - {£f {« 


Second Order Approximation, of Angular Velocity 


bl rCiT 
0) = {b> 


®1 + ®3®2 

02 - 636, 

83^ ®2®1 


{£'/ I [E] + 


0 

-e. 

0. 


0 

0 


= {bf {0} + 


0 0 

0 0 -0 

0 i 



{0} = {b}^([E] + [e^]){e} 


{0} 1 = {bf({0> + [e^Ko}) 


# pp • rp • 

[0a1 - 

BEAM LOCALLY ATTACHED TO INERTIAL COORDINATE TRANSEORMATION 
The bg axis is placed tangent to the centroid at each point. This 
accomplishes rotations 1 and 2. Rotation three is the torsional rota- 
tion, about the centroid (now b^) . Let T be the unit vector tangent 
to the centroid, located at (u^,U 2 »z) 


T = 


-1 ^2 -2 -3 

Vl + 


1A4 



Since T = ;b^, the direction cosines of must equal the components 
of T. 


= sin0„ = 


2 - “i/V^ 


C32 —sinO^cose^ = u* yWl + (up 2 + (up 2 


C33 = cos6^cos02 = 1 



1 + (up2 + (up 2 


The direction cosines may then be defined in terms of 


sin0- = - u' 
1 2 


' 2 

1 + (up^ 


COS0^ = 1 


'1 + (up^ 


sine 2 = 


u 


1 / + (“P" 


cosBj . 1 + (up^/ + (up^ + (up2 


The derivatives of 8^ and ©2 present in the angular velocity expression 
may be calculated by taking the derivatives sin0^ and sin02. The 
binomial theorem may then be used in order to evaluate the second-order 
approximation of the angular velocity. 


bi 


i-u: 




= {b}^ |[E] + 


1 

I t^u’ 

+ ' 

[ ^u' 

( 


0 



0 


0 


~u. 
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APPENDIX C 


VECTOR-DYADIC RELATIONSHIPS 

Dot and Cross Products 


U • V = {u}'^{v} 

■Q X V = [u] {v} 

where ja = 

V = 

{b> = 


/S fp 





[u] = 


0 -u. 


“^2 "l 


[u]{v> = - [v]{u} 

D "" u = ii>^[I][u]<b} 

where Q = {_b [ I ] 


[I] 


^11 ^12 




"13 


21 22 23 

J31 ^32 I33 


' ,;E BLAl^ NCff 
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APPENDIX D 


HYBRID COORDINATE tIETHOD USING ASSUMED 
MODE SHAPES FOR ELASTIC CONTINUA 

ABSTRACT ; The hybrid coordinate method provides equations of motion of 

1 

minimum dimension for a spacecraft with flexible appendages. Instead of 
the usual finite element approach, in which mode shapes are calculated from 
equations of vibration of the finite element assembly, this chapter provides 

' f 

an alternative formulation using assumed mode shapes. This proves useful 
for a class of simply modeled appendages for which mode shapes are provided 
by an outside agency, or are otherwise known. The results are shown to be 
compatible with the finite element formulation, as previously described. 
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I. INTRODUCTION 


The hybrid-coordinate method provides equations of motion of a space- 
craft with elastic (flexible) appendages. The appendages are modeled as an 
interconnected set of small rigid bodies interconnected by massless or mass- 
ive elastic bodies (finite elements). From the Newton-Euler approach, the 
equations of motion for each finite element and the rigid body portion of the 
spacecraft are formulated. The introduction of an appropriate coordinate 
transformation allows the finite element equations to be represented as 
decoupled vibration equations, which involve mode shapes and modal coordinates. . 
Since the vibration equations have been decoupled from each other, significant 
truncation of the higher order mode shapes can be accomplished. This leads 
to a set of equations where rotation of the rigid body portion of the space- 
craft is coupled to the vibration of the flexible appendages. These equations 
are of great practical use because the truncation procedure has significantly 
reduced the mnnber of degrees of freedom of the system without substantially 
sacrificing the fidelity of the results. 

The purpose of this chapter is to provide an alternative formulation 
for the hybrid coordinate method using assumed mode shapes. This approach 
will prove useful for simply modeled appendages. For these the mode shapes 
can be. determined from a continuum analysis using partial differential equa- 
tion methods. The truncation procedure is accomplished at the outset by 
eliminating the higher order modes of vibration. The equations of motion 
are formulated using a l^agrangian approach and the coordinate transformation 
is accomplished using the assumed mode shapes. The resulting equations of 
motion are then seen to be compatible with those arising from the finite 
element method. 

' blank- not 
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Figure 1. System Diagram. 

Center of mass of undeformed system (body fixed) 
Position of 0 at rest (inertially fixed) 

Center of mass of rigid body 

Center of mass of undeformed appendage 

Total System mass 

Appendage mass 

Appendage mass/length 

ConnectiCtn point of appendage 




II. MODEL 

The following derivation of equations of motion uses a model comprised 
of a central rigid body and a flexible cantilevered beam. Extensions to 
several appendages of arbitrary configurations may be made from the results 

of this simple model. The undeformed position of the appendage is taken to be 

constant relative to the rigid body. The transformation between the two is 
included in the derivation to facilitate the extension of the equations to 
cover a driven appendage. The angular rotations are assumed to be small as 

are the translational displacements. A diagram of the model is shown in Fig- 

ure 1. To summarize, the assuiuptions used in the following derivation are: 

• rigid body with cantilevered beam, 

• beam rest position constant relative to base, 

• small translations and rotations. 

No orthogonality requirements have been placed on the assumed mode 
shapes. The vibration equations are therefore coupled. Further coordinate 
transformations may be employed to decouple the vibration equations or to 
achieve vehicle normal modes, but the truncation procedure does not require 
this as it does with the finite element procedure. 

The vector bases employed in the derivation are: 

{i} : Inertially fixed basis 

{b} : Basis fixed in the rigid body 

{a} : Basis fixed in appendage prior to deformation 

where 

{b} - [BHi} 

[0] (E - 0 ) for small rotations 9 = 
la> - [clibl 

[c] *= constant for an undriven appendage. 
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For the undeformed system, the location of the center of mass is 
defined by 

/ p dm « 0 
SYS 

where £ is the generic position vector from the center of mass to the dif- 
ferential mass element. Evaluating this expression leads to 

^ L + m(r + r^) = 0 ^ 

where the quantities are shown in Figure 1 with L being the location of the 
rigid body center of mass and r^ the appendage center of mass. 

The dyadic of the undeformed system is defined by 

d” 

^AP-u 


n* r 

1° 

U “L 

RB 



^RB ^ 

-"rb 

n“=r 


U L 

AP-u 

AP* 


-a) - ■ (— a)(s+£a)] 

, III. THE LAGRANGIAN OF THE SYSTEM 

9 

The kinetic energy of the system is 

T =-|/*V'Vdm + V • V dm 
RB APP 

where V is the inertial velocity of a generic mass element. The kinetic 
energy for the rigid body yields 

I y’v • V dm - (^/-m)r®*R® + i w • . W 


where is the system mass, M .the appendage mass, and B is the center of mass 
of the rigid body with 
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?B = 5 - £ 


Rg *= Z - w X L. 


Expanding this expression and switching the reference point of the inertia 
dyadic to the system mass center gives 


i /* Y • V dm = ) Z«Z+ OD • n • to -L//-M)f*Z • (tOXL)] 

. RB RB ' 'L- -J 

The kinetic energy of the appendage is 

ly V . vam-iyi 


R • R dm 
-m -m 


where 


APP 


R^»Z + R + r + u 


•••00 / , 

+ H + S + E + R+r j . 

Vector differentiation with respect to the rotating reference frame is denoted 
by the "circle" above the vector. Here R and | are zero since they are fixed 
in the frame. The "dot" denotes differentiation relative to an inertial refer- 
ence frame. Expanding the expression and making use of the dyadic of the unde- 


appendage about the system mass centerYfl^ ) yields. 

'^AP-u' 

1 y'v • V dm = i MZ • z + 1 0 ) . [] 

1 • W + i /'u 

AP-u ^J~ 

APP 


* f u dm + Z • 

a X [m(r + ea)] 

APP 

Lx /-i 

+ ^ j^u • (0 X ^ R+r j 

dm 


udm 


APP 


Combining the terms for kinetic energy and eliminating terms produces 
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sits - • n 


^ J" u'u dm + Z • f u dm 
APP APP 

w»R x J" u dm + w • ^ (rxu) dm. 


The inertia dyadic of the undeformed system about the system center of mass is 
* 

lJ . The center of mass expression eliminated the term in the kinetic energy 
containing Z and w. 

From beam theory, the strain energy of the appendage is 


APP / 




The Lagrangian for the system is then 


*= Z*Z+-^(jO*Qj •o) + Y f u*u 

APP 

+ Z • u dm + 0) • R X J' u dm 


+ 03 • /* (rxu) dm - I f El (^\ • ( dr. 

/ \9r / 

APP APP ' / \ f 


The formation of the above has assumed that the undeformed appendage is 
fixed relative to the base ^ ~ j • next step to be taken is to 

assume small angle rotations and represent the Lagrangian in matrix form. The 
following matrices are used : 

Z »■ fi}^ {Z> 

03 “ (small rotations) 
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I 

I 

i; 

. K 

I 

‘i 

I 

u « {a}*^ {u} I 


g*- {b}^ I* {b} 
R = {b}^ {R> 
r = {b}^ {r} 



Retaining second order terms in the Lagrangian produces 

^ {Z} {0}^ I* {0} 

+ ^ J' {u}^ {ii} dm + {Z}[C] J* {u} dm 
APP . APP 

+ {0}"^ C ^ J' ^'^3' {0}^ C J* r {u} dm 

APP APP . 

t /* f( rp It 

- -2 ^ El {u } {u } dr. 

APP 
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Distributed coordinates are Introduced by the coordinate transformation 

n 

u(r,t) = 53 (t) 

1=1 ” 

where n Is the number of modes used to represent the displacement. In matrix 

form, the transformation Is 

n 

{u} = 53 

1=1 

■= [{>' I i ... i 

= [4>]{n} 

where [({)] is a 3 x n matrix with each column corresponding to a mode shape 
and {f)} contains n modal coordinates. 

This coordinate transformation yields the Lagrangian . 

5? = "I {Z}*^ {Z} + i {e}*^ I* {e} 

+ X2 {ri) + C {f|} 

+ {0}^ c(R x^+x^) ix]} - I {n}^ x^ {n> 

where 



APP 

APP 

X^ = jT r* [(|)] dm 
APP 

/ II ff II 

El [4> 3 [<t> 3 dr 

APP 
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The matrices and are of dimension 3 x n while the matrices X 2 and X^ 
are symmetric and of dimension nXn, The Lagrangian depends on n +6 generalized 
coordinates. Six coordinates describe the translation and rotation of the 
undeformed system and n modal coordinates describe the displacement from rest 
of the flexible appendage relative to the rigid base. 

IV. EQUATIONS OF MOTION 

The equations of motion for the system may now be derived from the 
Lagrangian in the traditional manner. The resulting n +6 equations may be 
represented in matrix form as 

at iz} + {t{} ^ 0 

I* {0} + (r x^+x^^Cii} = {t} 

X2 {ii} + {n> = {z} + ( xj R-X3){e} 

where {t} is the externally applied torque. The first matrix equation may be 
used to eliminate the translation from the vibration equations. This pro- 
duces n + 3 equations of the form 

I* {0} + ( R Xj^ + X 3 ){Ti> = {T} 

(*2 ’'l) - ^3 ) ® 

The matrices that provide coupling between the rotation and vibration 
in each equation may be seen to be transposes of each other. The equations 
can be written as 

• I* { 6 } - {fi} = {T} 

(*2 ^‘i) 

where 

fTi rsJ rp, 

6 X^ R - X 3 . 
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V. COMPATIBILITY WITH FINITE 


ELEMENT EQUATIONS OF MOTION 


The equations of motion derived from the continuum ana y sis are similar 
in structure to those derived from the finite element analysis (Ref. 1). The 
differences between the two appear in the assumptions made in the continuum 
analysis: 

• No orthogonality properties 
’• Rotary inertia effect ignored ^ 

(No differential rotation of appendage mass 
elements due to deformation) 

Orthogonality properties can be applied to the continuum analysis 
vibration equations by a suitable coordinate transformation. The orthogonality 
properties are not needed to permit truncation as is the case in the finite 
element analysis. 


The equations of motion from a finite element analysis are shown by 
equations (287) to (289) of Reference 1, 

*VirT~ 

I 6 - 5 n = T * 

•• —.1 —2— — 

p + 2 ^ ot \ + an = 66 


6 = -(j) 


m(e - £ r-7£ ) 

\ OE EO EO/ 


The overbar indicates truncation. If the damping is eliminated and 
the orthogonality condition relaxed (after truncation), the equations become 

^ — np 

I 8 - 6 n = T 


where 


Ti I* T t - — 

M (p n + (p K (p n = <5 d 

M = M ^ E - iCgQ 

f 

K = Stiffnes- matrix. 
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REPRODUCIBILITY OF THE 
ORK^N^AL PAGE IS POOB 



Normally, the coordinate transformation <J> includes mode shapes with 
translation and rotation of the finite elements. To agree with the continuum 
analysis, no rotations of the finite elements will be allowed. The coordinate 
transformation <P will then be a 6n x N matrix represented by 


where 


5 - 




{<!}’ 


> 

n 


With the above limitations, the matrix multiplication can be performed 
in the finite element equations and the terms may be compared with those 
from the continuum analysis. 

For the augmented mass matrix, the finite element analysis results in 

n 




4 

T 

4 

k 


k 


n ' n 

+ E S 


k l 

m m 


Z=1 k=l I 


D, 


This is compatible with the result from the continuum analysis 

ij APP L J L J 


// 

APP ■ 


m(r)m(r ) 




dr dr. 
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For the 6 matrix, the finite element analysis results in 


[J], 


m. 

x=l X 




(R . r,) 


This is compatible with the result from the continuum analysis 

T 


[<5]j m(r) ^ R + r j dr 


Thus, if the number of finite elements were increased without limit, the 
finite element equations would be identical to the continuum analysis 
equations. 


VI. CONCLUSION 


With the foregoing results it becomes possible to accomplish a hybrid 
coordinate dynamic analysis for a system with appendages defined only in 
terms of modal data based on a continuum analysis. 


I 
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